rm(list=ls())
require(stringr)
require(zoo)
require(openair)
source("concord_filter.R")
source("cospectra_plot3.R")
na.count<-function(x) sum(is.na(x))
na.mean<-function(x) ifelse(is.nan(mean(x,na.rm=T)),NA,mean(x,na.rm=T))
library(REddyProc)
library(car)
library(tidyverse)
library(ggpubr)
library(rstatix)
library(broom)
library(psych)
library(lme4)
library(modelsummary)
hc: define all I/O upfront and reuse them for consistency
## change root.path as needed
root.path<-"C:\\Users\\tfens\\R_REPOS\\Flux_processing\\Concord_R_Code\\Concord_Post_Process\\"
#root.path<-"D:\\Housen\\Flux\\Data-exploring\\02_Concord\\"
eddypro.path<-paste0(root.path,"01_data\\Master_MET_Eddy\\")## this is where the EddyPro outputs located
path.in_met<-paste0(root.path,"01_data\\Master_MET_Eddy\\")
plot.path<-paste0(root.path,"02_output\\02_spectrum_plot\\") ## this is where the output plots located
# hc: Use this for storing combined file
path.out<-paste0(root.path,"02_output\\03_Met_Eddy_combined\\",sep="")
hc: Keep anything that needs to change regularly here
## use follows to specify the versions of EddyPro outputs
# file name of the master EddyPro file
cdata.file2<-paste0("master_eddy_pro_concord.csv")
## Use follows to define binned spectrum file names
cdata.proc.time<-"2020-04-17T161020"
cdata.proc.ext<-"_adv"
#20190614-1700_binned_cospectra_2020-04-17T161020_adv
## use follows to specify the versions of master met file
file.name<-paste0("MET_data_master.csv")
#"C:\Users\Tommy\flux\Data-exploring\02_Concord\01_Proccessed_Data\Master_Eddy"
#"C:\Users\Tommy\flux\Data-exploring\02_Concord\01_Proccessed_Data\Master_Eddy\master_eddy_pro_concord.csv"
#"C:\Users\Tommy\flux\Data-exploring\02_Concord\01_Proccessed_Data\Master_Eddy\eddypro_binned_cospectra"
#"C:\Users\Tommy\flux\Data-exploring\02_Concord\01_Proccessed_Data\Master_Eddy\eddypro_full_cospectra"
#"C:\Users\Tommy\flux\Data-exploring\02_Concord\01_Proccessed_Data\Master_Eddy\spectrum_plot"
#parse variable names and define N/As #remove time periods that does not have enough record of high frequency data
## read in full output file
cdata<-read.csv(paste(eddypro.path,cdata.file2,sep=""),
skip=3,
header=F,
na.strings="-9999",
stringsAsFactors=F)
colnames(cdata)<-colnames(read.csv(paste(eddypro.path,cdata.file2,sep=""),
skip=1,
header=T,
na.strings="NaN"))
cdata<-cdata[cdata$filename!="not_enough_data"&
!is.na(cdata$used_records),]
head(cdata)
tail(cdata)
NA
met_data_master<-read.csv(paste(path.in_met,file.name,sep=""),
header=F,
skip=4,
na.strings=c("NAN","-7999"),
stringsAsFactors = F)
colnames(met_data_master)<-colnames(
read.csv(paste(path.in_met,file.name,sep=""),
header=T,
skip= 1))
head(met_data_master)
tail(met_data_master)
NA
NA
NA
NA
cdata$TIMESTAMP<-strptime(paste(cdata$date,cdata$time,sep=" "),
format="%m/%d/%Y %H:%M",
tz = "Etc/GMT-8")
#cdata$TIMESTAMP=cdata$TIMESTAMP+1800 don't need to add 1800. already in endtime format
#cdata$TIMESTAMP=cdata$TIMESTAMP #just rename it end. so that we know
cdata$time.id<-cdata$TIMESTAMP$year+1900+
(cdata$TIMESTAMP$yday)/366+
(cdata$TIMESTAMP$hour)/366/24+
(cdata$TIMESTAMP$min)/366/24/60
cdata$time.id[1:50]
[1] 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495
[13] 2019.495 2019.495 2019.495 2019.495 2019.495 2019.496 2019.496 2019.496 2019.496 2019.496 2019.496 2019.496
[25] 2019.496 2019.496 2019.496 2019.496 2019.496 2019.496 2019.496 2019.496 2019.496 2019.496 2019.496 2019.497
[37] 2019.497 2019.497 2019.497 2019.497 2019.497 2019.497 2019.497 2019.497 2019.497 2019.497 2019.497 2019.497
[49] 2019.497 2019.497
plot(cdata$TIMESTAMP,cdata$time.id)
which(duplicated(cdata$time.id))
integer(0)
Taking the met_data and turning the time stamp into posixt format creating a time id for the MET Data so I I can join the MET and Eddy Pro Data
met_data_master$TIMESTAMP<-strptime(met_data_master$TIMESTAMP,
format ="%m/%d/%Y %H:%M",
tz = "Etc/GMT-8")
met_data_master$time.id <-met_data_master$TIMESTAMP$year+1900+
(met_data_master$TIMESTAMP$yday)/366+
(met_data_master$TIMESTAMP$hour)/366/24+
(met_data_master$TIMESTAMP$min)/366/24/60
met_data_master$time.id[1:20]
[1] 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495 2019.495
[13] 2019.495 2019.495 2019.495 2019.495 2019.495 2019.496 2019.496 2019.496
plot(met_data_master$TIMESTAMP,met_data_master$time.id)
which(duplicated(met_data_master$time.id))
[1] 42213 42214 42215 42216 42217 42218 42219 42220 42221 42222 42223 42224 42225 42226 42227 42228 42229 42230 42231
[20] 42232 42233 42234 42235 42236 42237
#removing duplicated rows
met_data_master<-met_data_master[!duplicated(met_data_master$time.id), ]
#checking that this worked
which(duplicated(met_data_master$time.id))
integer(0)
head(met_data_master)
met_data_master$TIMESTAMP[1:20]
[1] "2019-07-01 00:00:00 +08" "2019-07-01 00:30:00 +08" "2019-07-01 01:00:00 +08" "2019-07-01 01:30:00 +08"
[5] "2019-07-01 02:00:00 +08" "2019-07-01 02:30:00 +08" "2019-07-01 03:00:00 +08" "2019-07-01 03:30:00 +08"
[9] "2019-07-01 04:00:00 +08" "2019-07-01 04:30:00 +08" "2019-07-01 05:00:00 +08" "2019-07-01 05:30:00 +08"
[13] "2019-07-01 06:00:00 +08" "2019-07-01 06:30:00 +08" "2019-07-01 07:00:00 +08" "2019-07-01 07:30:00 +08"
[17] "2019-07-01 08:00:00 +08" "2019-07-01 08:30:00 +08" "2019-07-01 09:00:00 +08" "2019-07-01 09:30:00 +08"
Create a full timestamp based on the earliest and latest timestamps in EddyPro and Met files Use this full timestamp later when merging EddyPro and master met files
# create a full timestamp, 30 mins
full.time<-data.frame(TIMESTAMP=
seq.POSIXt(min(min(met_data_master$TIMESTAMP),min(cdata$TIMESTAMP)),
max(max(met_data_master$TIMESTAMP),max(cdata$TIMESTAMP)),
units = "seconds", by = 1800),
stringsAsFactors=F)
full.time$TIMESTAMP<-strptime(full.time$TIMESTAMP,
format ="%Y-%m-%d %H:%M:%S",
tz = "Etc/GMT-8")
full.time$time.id <-full.time$TIMESTAMP$year+1900+
(full.time$TIMESTAMP$yday)/366+
(full.time$TIMESTAMP$hour)/366/24+
(full.time$TIMESTAMP$min)/366/24/60
print(paste("Starting timestamp:",full.time$TIMESTAMP[1]))
[1] "Starting timestamp: 2019-07-01"
print(paste("Ending timestamp:",full.time$TIMESTAMP[nrow(full.time)]))
[1] "Ending timestamp: 2022-01-28 10:30:00"
head(full.time)
NA
NA
#Joining the Met_Data and Eddy Pro Data Sets using time stamp from the full.time dataframe Also create a doy.id, unique for each date, later used in aggregating daily values
cdata<- merge.data.frame(full.time,
cdata[,-which(colnames(cdata)=="TIMESTAMP")],
by = "time.id",
all = TRUE,
sort = TRUE)
#all=true what ever appears in each file is show in the file data file. sort tries to sort each data frame by merging. in this case probably doesnot matter. timestamp
cdata<- merge.data.frame(cdata,
met_data_master[,-which(colnames(met_data_master)=="TIMESTAMP")],
by = "time.id",
all = TRUE,
sort = TRUE)
#seeing which rows are duplicated
which(duplicated(cdata$TIMESTAMP))
integer(0)
#removing duplicated rows
# cdata<-cdata[!duplicated(cdata$TIMESTAMP), ]
#
# #checking that this worked
# which(duplicated(cdata$TIMESTAMP))
## create a unique DOY id
cdata$doy.id<-full.time$TIMESTAMP$year+1900+
(full.time$TIMESTAMP$yday)/366
colnames(cdata)
[1] "time.id" "TIMESTAMP" "filename"
[4] "date" "time" "DOY"
[7] "daytime" "file_records" "used_records"
[10] "Tau" "qc_Tau" "H"
[13] "qc_H" "LE" "qc_LE"
[16] "co2_flux" "qc_co2_flux" "h2o_flux"
[19] "qc_h2o_flux" "H_strg" "LE_strg"
[22] "co2_strg" "h2o_strg" "co2_v.adv"
[25] "h2o_v.adv" "co2_molar_density" "co2_mole_fraction"
[28] "co2_mixing_ratio" "co2_time_lag" "co2_def_timelag"
[31] "h2o_molar_density" "h2o_mole_fraction" "h2o_mixing_ratio"
[34] "h2o_time_lag" "h2o_def_timelag" "sonic_temperature"
[37] "air_temperature" "air_pressure" "air_density"
[40] "air_heat_capacity" "air_molar_volume" "ET"
[43] "water_vapor_density" "e" "es"
[46] "specific_humidity" "RH" "VPD"
[49] "Tdew" "u_unrot" "v_unrot"
[52] "w_unrot" "u_rot" "v_rot"
[55] "w_rot" "wind_speed" "max_wind_speed"
[58] "wind_dir" "yaw" "pitch"
[61] "roll" "u." "TKE"
[64] "L" "X.z.d..L" "bowen_ratio"
[67] "T." "model" "x_peak"
[70] "x_offset" "x_10." "x_30."
[73] "x_50." "x_70." "x_90."
[76] "un_Tau" "Tau_scf" "un_H"
[79] "H_scf" "un_LE" "LE_scf"
[82] "un_co2_flux" "co2_scf" "un_h2o_flux"
[85] "h2o_scf" "spikes_hf" "amplitude_resolution_hf"
[88] "drop_out_hf" "absolute_limits_hf" "skewness_kurtosis_hf"
[91] "skewness_kurtosis_sf" "discontinuities_hf" "discontinuities_sf"
[94] "timelag_hf" "timelag_sf" "attack_angle_hf"
[97] "non_steady_wind_hf" "u_spikes" "v_spikes"
[100] "w_spikes" "ts_spikes" "co2_spikes"
[103] "h2o_spikes" "chopper_LI.7500" "detector_LI.7500"
[106] "pll_LI.7500" "sync_LI.7500" "mean_value_RSSI_LI.7500"
[109] "u_var" "v_var" "w_var"
[112] "ts_var" "co2_var" "h2o_var"
[115] "w.ts_cov" "w.co2_cov" "w.h2o_cov"
[118] "vin_sf_mean" "co2_mean" "h2o_mean"
[121] "dew_point_mean" "co2_signal_strength_7500_mean" "RECORD"
[124] "BattV_Avg" "PTemp_C_Avg" "AM25T_ref_Avg"
[127] "TC_Avg.1." "TC_Avg.2." "TC_Avg.3."
[130] "TC_Avg.4." "TC_Avg.5." "TC_Avg.6."
[133] "TC_Avg.7." "TC_Avg.8." "TC_Avg.9."
[136] "TC_Avg.10." "TC_Avg.11." "TC_Avg.12."
[139] "TC_Avg.13." "TC_Avg.14." "TC_Avg.15."
[142] "TC_Avg.16._control" "TC_Avg.17._control" "TC_Avg.18._control"
[145] "TC_Avg.19._control" "TC_Avg.20._control" "TC_Avg.21."
[148] "TC_Avg.22." "TC_Avg.23." "TC_Avg.24."
[151] "TC_Avg.25." "AirT_Avg" "RH_Avg"
[154] "AtmPressure_Avg" "NR_mV_Avg" "NR_Wm2_Avg"
[157] "PAR_in_mV_Avg" "PAR_in_uEm2_Avg" "PAR_out_mV_Avg"
[160] "PAR_out_uEm2_Avg" "SHF_1_mV_Avg" "SHF_1_Wm2_Avg"
[163] "SHF_2_mV_Avg" "SHF_2_Wm2_Avg" "WaterP_Avg"
[166] "WaterT_Avg" "Precip_mm_Tot" "VWC_Avg"
[169] "EC_Avg" "T_Avg" "P_Avg"
[172] "PA_Avg" "VR_Avg" "VWC_2_Avg_Control"
[175] "EC_2_Avg_Control" "T_2_Avg_Control" "P_2_Avg_Control"
[178] "PA_2_Avg_Control" "VR_2_Avg_Control" "doy.id"
cdata
NA
#Creating a CSV File of my combined Master File!# #hc: modify the filename starting with Date when the file is created
write.csv(cdata,
paste(path.out,Sys.Date(),"_master_eddy_met_concord_prefiltering.csv",sep=""),
quote = T,
row.names = F)
#Filter data based on predefined criteria #filtering criteria are defined in concord_filter function # A few variables also convert Unit
# plot(
# cdata$Correct_NR[cdata$wind_dir >= 230 &
# cdata$wind_dir <= 300 &
# cdata$daytime == 0 &
# cdata$u. > 0.1],
# cdata$co2_flux[cdata$wind_dir >= 230 &
# cdata$wind_dir <= 300 & cdata$daytime == 0 & cdata$u. > 0.1],
# # xaxt= 'n',
# # yaxt= 'n',
# abline(h = 0, col = "darkgrey"),
# ylim = c(-20, 20),
# xlim = c(0, 600),
# xlab = 'Net Radiation',
# ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
# -2
# } ~ s ^ {
# -1
# } ~ ')'),
#
# #main='Night time fluxes by Net Radiation',
#
# cex = 0.6,
# col = "red"
# )
# plot(cdata$RH)
#
# plot(cdata$RH_Avg)
#
# plot(cdata$VPD)
#
# plot(cdata$time.id, cdata$PAR_in_w_m_2
# ylim = c(0,750),
# xlim = c(2019.928, 2020.9)
#)
cdata<-concord_filter(data.in=cdata)
colnames(cdata)
[1] "time.id" "TIMESTAMP" "filename"
[4] "date" "time" "DOY"
[7] "daytime" "file_records" "used_records"
[10] "Tau" "qc_Tau" "H"
[13] "qc_H" "LE" "qc_LE"
[16] "co2_flux" "qc_co2_flux" "h2o_flux"
[19] "qc_h2o_flux" "H_strg" "LE_strg"
[22] "co2_strg" "h2o_strg" "co2_v.adv"
[25] "h2o_v.adv" "co2_molar_density" "co2_mole_fraction"
[28] "co2_mixing_ratio" "co2_time_lag" "co2_def_timelag"
[31] "h2o_molar_density" "h2o_mole_fraction" "h2o_mixing_ratio"
[34] "h2o_time_lag" "h2o_def_timelag" "sonic_temperature"
[37] "air_temperature" "air_pressure" "air_density"
[40] "air_heat_capacity" "air_molar_volume" "ET"
[43] "water_vapor_density" "e" "es"
[46] "specific_humidity" "RH" "VPD"
[49] "Tdew" "u_unrot" "v_unrot"
[52] "w_unrot" "u_rot" "v_rot"
[55] "w_rot" "wind_speed" "max_wind_speed"
[58] "wind_dir" "yaw" "pitch"
[61] "roll" "u." "TKE"
[64] "L" "X.z.d..L" "bowen_ratio"
[67] "T." "model" "x_peak"
[70] "x_offset" "x_10." "x_30."
[73] "x_50." "x_70." "x_90."
[76] "un_Tau" "Tau_scf" "un_H"
[79] "H_scf" "un_LE" "LE_scf"
[82] "un_co2_flux" "co2_scf" "un_h2o_flux"
[85] "h2o_scf" "spikes_hf" "amplitude_resolution_hf"
[88] "drop_out_hf" "absolute_limits_hf" "skewness_kurtosis_hf"
[91] "skewness_kurtosis_sf" "discontinuities_hf" "discontinuities_sf"
[94] "timelag_hf" "timelag_sf" "attack_angle_hf"
[97] "non_steady_wind_hf" "u_spikes" "v_spikes"
[100] "w_spikes" "ts_spikes" "co2_spikes"
[103] "h2o_spikes" "chopper_LI.7500" "detector_LI.7500"
[106] "pll_LI.7500" "sync_LI.7500" "mean_value_RSSI_LI.7500"
[109] "u_var" "v_var" "w_var"
[112] "ts_var" "co2_var" "h2o_var"
[115] "w.ts_cov" "w.co2_cov" "w.h2o_cov"
[118] "vin_sf_mean" "co2_mean" "h2o_mean"
[121] "dew_point_mean" "co2_signal_strength_7500_mean" "RECORD"
[124] "BattV_Avg" "PTemp_C_Avg" "AM25T_ref_Avg"
[127] "TC_Avg.1." "TC_Avg.2." "TC_Avg.3."
[130] "TC_Avg.4." "TC_Avg.5." "TC_Avg.6."
[133] "TC_Avg.7." "TC_Avg.8." "TC_Avg.9."
[136] "TC_Avg.10." "TC_Avg.11." "TC_Avg.12."
[139] "TC_Avg.13." "TC_Avg.14." "TC_Avg.15."
[142] "TC_Avg.16._control" "TC_Avg.17._control" "TC_Avg.18._control"
[145] "TC_Avg.19._control" "TC_Avg.20._control" "TC_Avg.21."
[148] "TC_Avg.22." "TC_Avg.23." "TC_Avg.24."
[151] "TC_Avg.25." "AirT_Avg" "RH_Avg"
[154] "AtmPressure_Avg" "NR_mV_Avg" "NR_Wm2_Avg"
[157] "PAR_in_mV_Avg" "PAR_in_uEm2_Avg" "PAR_out_mV_Avg"
[160] "PAR_out_uEm2_Avg" "SHF_1_mV_Avg" "SHF_1_Wm2_Avg"
[163] "SHF_2_mV_Avg" "SHF_2_Wm2_Avg" "WaterP_Avg"
[166] "WaterT_Avg" "Precip_mm_Tot" "VWC_Avg"
[169] "EC_Avg" "T_Avg" "P_Avg"
[172] "PA_Avg" "VR_Avg" "VWC_2_Avg_Control"
[175] "EC_2_Avg_Control" "T_2_Avg_Control" "P_2_Avg_Control"
[178] "PA_2_Avg_Control" "VR_2_Avg_Control" "doy.id"
[181] "air_temperature_adj" "Correct_NR" "Rg"
[184] "out_going_rad" "Correct_shf_1" "Correct_shf_2"
#Creating a CSV File of combined Master File (post-filtering)
write.csv(cdata,
paste(path.out,Sys.Date(),"_master_eddy_met_concord_postfiltering.csv",sep=""),
quote = T,
row.names = F)
#Generate simple timeseries plots post-filtered time series plot per variable LOESS fit and daily average are plotted to show temporal dynamics
#summary(cdata)
##############################################################
## Generic time series plot
target.plot.var<-c("co2_flux","LE","H","u.",
"co2_mixing_ratio","h2o_mixing_ratio",
"air_temperature_adj","RH_Avg",
"wind_speed","wind_dir","mean_value_RSSI_LI.7500",
"AirT_Avg", "Correct_NR", "Correct_shf_1" , "Correct_shf_2",
"VWC_Avg","Precip_mm_Tot" ,
"PAR_in_uEm2_Avg","PAR_out_uEm2_Avg",
"Rg",
"AtmPressure_Avg",
"TC_Avg.1.","TC_Avg.2.","TC_Avg.3.","TC_Avg.4.","TC_Avg.5.")
target.plot.var.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(LE~'('~W~m^{-2}~')'),
expression(H~'('~W~m^{-2}~')'),
expression(u['*']~'('~m~s^{-1}~')'),
expression(CO[2]~'('~ppm~')'),
expression(Water~vapor~'('~ppt~')'),
expression(Air~temperature~'('~degree~C~')'),
expression(Relative~humidity~'('~percent~')'),
expression(Wind~speed~'('~m~s^{-1}~')'),
expression(Wind~direction~'('~degree~')'),
expression(LI7500~signal~strengh ~'('~'-'~')'),
expression(Air~temperature~MET~'('~degree~C~')'),
expression(Net~Radiation~'('~W~m^{-2}~')'),
expression(Soil~Heatflux1~'('~W~m^{-2}~')'),
expression(Soil~Heatflux2~'('~W~m^{-2}~')'),
expression(Soil~water~content~'('~m^{3}~m^{-3}~')'),
expression(Precipitation~'('~mm~')'),
expression(Incoming~PAR~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(Outgoing~PAR~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(Rg~'('~W~m^{-2}~')'),
expression(Atmospheric~pressure~'('~kPa~')'),
expression(Soil~temperature~1~'('~degree~C~')'),
expression(Soil~temperature~2~'('~degree~C~')'),
expression(Soil~temperature~3~'('~degree~C~')'),
expression(Soil~temperature~4~'('~degree~C~')'),
expression(Soil~temperature~5~'('~degree~C~')'))
for(k1 in 1:length(target.plot.var)){
## locate the start of each month
month.loc<-which(cdata$TIMESTAMP$mday==1&
cdata$TIMESTAMP$hour==0&
cdata$TIMESTAMP$min==0)
month.ticks <- seq(cdata$TIMESTAMP[month.loc[1]],
cdata$TIMESTAMP[month.loc[length(month.loc)]],by="months")
png(paste0(plot.path,
"Concord_",
cdata$TIMESTAMP$year[1]+1900,"_",
cdata$TIMESTAMP$yday[1]+1,"_",
cdata$TIMESTAMP$year[nrow(cdata)]+1900,"_",
cdata$TIMESTAMP$yday[nrow(cdata)]+1,"_",
target.plot.var[k1],"_",
Sys.Date(),".png"),
width=6,
height=4,
units="in",
res=300,
pointsize = 11,
bg = "white")
par(oma=c(0.5,0.5,0.5,0.5),mar=c(4,4.5,0,0))
plot(cdata$TIMESTAMP,
cdata[,target.plot.var[k1]],
xlab="TIMESTAMP",
ylab=target.plot.var.title[k1],
cex=0.5,col="grey",bg="lightgrey",
las=1,pch=21,
xaxs="i",yaxs="i"
)
abline(h=0,col="darkgrey")
## daily averge line
daily.tmp<-data.frame(date=tapply(cdata$time.id,cdata$doy.id,min),
daily=tapply(cdata[,target.plot.var[k1]],cdata$doy.id,na.mean))
lines(cdata$TIMESTAMP[which(cdata$time.id %in% daily.tmp$date)],
daily.tmp$daily,
lwd=1.5,col="black")
## loess fit line
loess.tmp<-loess(cdata[,target.plot.var[k1]]~c(1:nrow(cdata)),span=0.1)
## avoid plotting long gaps
long.gap<-which(loess.tmp$x[-1]-loess.tmp$x[-length(loess.tmp$x)]>48*3)
if(length(long.gap)>0){ ## skip if any long gaps
if(long.gap[1]>1) long.gap<-c(1,long.gap)
if(long.gap[length(long.gap)]<length(loess.tmp$x)) long.gap<-c(long.gap,length(loess.tmp$x))
for(k1 in 1:(length(long.gap)-1)){
lines(cdata$TIMESTAMP[round(loess.tmp$x[c((long.gap[k1]+1):(long.gap[k1+1]-1))])],
loess.tmp$fitted[c((long.gap[k1]+1):(long.gap[k1+1]-1))],
lwd=1.5,col="red")
}
}else{
lines(cdata$TIMESTAMP[round(loess.tmp$x)],
loess.tmp$fitted,
lwd=1.5,col="red")
}
axis(1, at = month.ticks, labels = FALSE, tcl = -0.3)
dev.off()
}
#Generate timeseries plots color-coded by wind direction similar plot as previous time series plot Color-coded data points based on wind direction, currently 2 groups hc: Revise thresholds for WD groups if needed
##############################################################
## Generic time series plot
WD.grp<-rep(2,nrow(cdata))
WD.grp[which(!is.na(cdata$wind_dir)&
(cdata$wind_dir>120&cdata$wind_dir<=300))]<-1
WD.legend<-c("Southwest wind","Northeast wind")
WD.col<-c(rgb(0,0,1,0.5,maxColorValue=1),rgb(1,0,0,0.5,maxColorValue=1))
target.plot.var<-c("co2_flux","LE","H","u.")
target.plot.var.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(LE~'('~W~m^{-2}~')'),
expression(H~'('~W~m^{-2}~')'),
expression(u['*']~'('~m~s^{-1}~')'))
for(k1 in 1:length(target.plot.var)){
## locate the start of each month
month.loc<-which(cdata$TIMESTAMP$mday==1&
cdata$TIMESTAMP$hour==0&
cdata$TIMESTAMP$min==0)
month.ticks <- seq(cdata$TIMESTAMP[month.loc[1]],
cdata$TIMESTAMP[month.loc[length(month.loc)]],by="months")
png(paste0(plot.path,
"Concord_",
cdata$TIMESTAMP$year[1]+1900,"_",
cdata$TIMESTAMP$yday[1]+1,"_",
cdata$TIMESTAMP$year[nrow(cdata)]+1900,"_",
cdata$TIMESTAMP$yday[nrow(cdata)]+1,"_",
target.plot.var[k1],"_color_",
Sys.Date(),".png"),
width=8,
height=4,
units="in",
res=300,
pointsize = 11,
bg = "white")
par(oma=c(4,4.5,0.5,0.5),mar=c(0,0.5,0,0.5),fig=c(0,0.7,0,1))
plot(cdata$TIMESTAMP,
cdata[,target.plot.var[k1]],
xlab="",
ylab="",
cex=0.5,col=WD.col[WD.grp],
las=1,pch=16,
xaxs="i",yaxs="i"
)
mtext(side=2,target.plot.var.title[k1],line=3)
mtext(side=1,"TIMESTAMP",line=2.8)
abline(h=0,col="darkgrey")
axis(1, at = month.ticks, labels = FALSE, tcl = -0.3)
par(fig=c(0.7,1,0,1),new=T)
hist0<-hist(cdata[,target.plot.var[k1]],
plot=F,nclass=50)
hist1<-hist(cdata[,target.plot.var[k1]][WD.grp==1],
plot=F,breaks=hist0$breaks)
hist2<-hist(cdata[,target.plot.var[k1]][WD.grp==2],
plot=F,breaks=hist0$breaks)
barplot(hist1$counts,
axes=F,
horiz=T,
ylim=c(0,length(hist1$breaks)+1),
xlim=c(-5,max(c(hist1$counts,hist2$counts))),
space=0,col=WD.col[1],border=NA) # barplot
barplot(hist2$counts,
axes=F,
add=T,
horiz=T,
ylim=c(0,length(hist1$breaks)+1),
xlim=c(-5,max(c(hist1$counts,hist2$counts))),
space=0,col=WD.col[2],border=NA) # barplot
legend(0,
length(hist1$breaks)+1,
fill=WD.col,border=NA,
legend=WD.legend,bty="n",
cex=0.9)
dev.off()
}
#Windrose plot A simple windrose plot, using openair package
png(paste0(plot.path,
"Concord_",
cdata$TIMESTAMP$year[1]+1900,"_",
cdata$TIMESTAMP$yday[1]+1,"_",
cdata$TIMESTAMP$year[nrow(cdata)]+1900,"_",
cdata$TIMESTAMP$yday[nrow(cdata)]+1,"_",
"Windrose_",
Sys.Date(),".png"),
width=5,height=5,units="in",res=300,pointsize=10)
openair::windRose(mydata=cdata[,c("wind_speed","wind_dir")],
ws="wind_speed",
wd="wind_dir",
ws.int = 0.5, angle = 15,
dig.lab=2,
paddle=F,key.position = "bottom")
dev.off()
null device
1
```r
colnames(cdata)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#Net Exchange of Co2 by U*. filter by night time and day time
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjI1xuI0RheXRpbWUgYW5kIG5pZ2h0dGltZSBwbG90XG5zY2F0dGVyLnNtb290aCggY2RhdGEkdS4sIGNkYXRhJGNvMl9mbHV4LFxuICAgICAgICAgICAgICAgeWxhYj1leHByZXNzaW9uKEZDficoJ35tdX5tb2x+bV57LTJ9fnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIHhsYWI9IGV4cHJlc3Npb24odVsnKiddficoJ35tfnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIG1haW49J0RheXRpbWUgYW5kIE5pZ2h0dGltZScpXG5cbiNkYXkgYW5kIG5pZ2h0IGxpbmVhciByZWxhdGlvbnNoaXBcbmxtX2FsbF90aW1lPC0gbG0odS4gflxuICAgICAgICAgICAgICAgICAgIGNvMl9mbHV4ICBcbiAgICAgICAgICAgICAgICAgICBcbiwgXG5kYXRhID0gY2RhdGEpXG5cbnN1bW1hcnkobG1fYWxsX3RpbWUpXG5cbiMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjI1xuI0RheXRpbWUgcGxvdFxuc2NhdHRlci5zbW9vdGgoc3Vic2V0KGNkYXRhLCBkYXl0aW1lPT0nMScpJHUuLFxuICAgICBzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScxJykkY28yX2ZsdXgsIFxuICAgICBcbiAgICAgICAgICAgICAgIHlsYWI9ZXhwcmVzc2lvbihGQ34nKCd+bXV+bW9sfm1eey0yfX5zXnstMX1+JyknKSwgXG4gICAgICAgICAgICAgICB4bGFiPSBleHByZXNzaW9uKHVbJyonXX4nKCd+bX5zXnstMX1+JyknKSwgXG4gICAgICAgICAgICAgICBtYWluPSdEYXl0aW1lICcpXG5cbiMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjI1xuI2RheXRpbWUgbGluZWFyIHJlbGF0aW9uc2hpcCBiZXR3ZWVuIGNvMiBmbHV4IGFuZCBVJlxubG1fZGF5X3RpbWU8LSBsbShzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScxJykkdS4gflxuICAgICAgICAgICAgICAgICAgIHN1YnNldChjZGF0YSwgZGF5dGltZT09JzEnKSRjbzJfZmx1eCBcbiAgICAgICAgICAgICAgICAgICBcbiwgXG5kYXRhID0gY2RhdGEpXG5cbnN1bW1hcnkobG1fZGF5X3RpbWUpXG4jIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjI1xuI05pZ2h0dGltZSBwbG90IGZlYl9tYXJjaFxuc21vb3RoU2NhdHRlcihzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgRE9ZPjMxICYgRE9ZIDwxMjApJHUuLFxuICAgICBzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgRE9ZPjMxICYgRE9ZIDwxMjAgKSRjbzJfZmx1eCwgXG4gICAgIFxuICAgICAgICAgICAgICAgeWxhYj1leHByZXNzaW9uKEZDficoJ35tdX5tb2x+bV57LTJ9fnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIHhsYWI9IGV4cHJlc3Npb24odVsnKiddficoJ35tfnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIG1haW49J05pZ2h0dGltZSAnLFxuICAgICAgY2V4PSAwLjUpXG5cbiNuaWdodHRpbWUgbGluZWFyIHJlbGF0aW9uc2hpcCBiZXR3ZWVuIGNvMiBmbHV4IGFuZCBVJlxubG1fbmlnaHRfdGltZV9mZWJfbWFyPC0gbG0oc3Vic2V0KGNkYXRhLCBkYXl0aW1lPT0nMCcmIERPWT4zMSAmIERPWSA8MTIwKSR1Ln5cbiAgICAgc3Vic2V0KGNkYXRhLCBkYXl0aW1lPT0nMCcmIERPWT4zMSAmIERPWSA8MTIwICkkY28yX2ZsdXhcbiwgZGF0YSA9IGNkYXRhKVxuXG5zdW1tYXJ5KGxtX25pZ2h0X3RpbWVfZmViX21hcilcblxuXG5cbiMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjI1xuI05pZ2h0dGltZSBwbG90IERyeSBTZWFzb24uIE1heSAxKDEyMSkgLU5vdiAxKDMwNSlcbnNjYXR0ZXIuc21vb3RoKHN1YnNldChjZGF0YSwgZGF5dGltZT09JzAnJiBET1k+MTIwICYgRE9ZIDwzMDYpJHUuLFxuICAgICBzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgRE9ZPjEyMCAmIERPWSA8MzA2ICkkY28yX2ZsdXgsIFxuICAgICBcbiAgICAgICAgICAgICAgIHlsYWI9ZXhwcmVzc2lvbihGQ34nKCd+bXV+bW9sfm1eey0yfX5zXnstMX1+JyknKSwgXG4gICAgICAgICAgICAgICB4bGFiPSBleHByZXNzaW9uKHVbJyonXX4nKCd+bX5zXnstMX1+JyknKSxcbiAgICAgeWxpbSA9IGMoMCwxMCksXG4gICAgICAgICAgICAgICBtYWluPSdOaWdodHRpbWUgRHJ5IFNlYXNvbiAoTWF5IDEtIE5vdi4gMSkgJyxcbiAgICAgIGNleD0gMC41KVxuXG4jIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyNcbiNOaWdodHRpbWUgcGxvdCBXZXQgU2Vhc29uLCBOb3YgMi0gQXByaWwgMzBcblxuc2NhdHRlci5zbW9vdGgoc3Vic2V0KGNkYXRhLCBkYXl0aW1lPT0nMCcmIERPWT4zMDUgfCBET1kgPDEyMSkkdS4sXG4gICAgIHN1YnNldChjZGF0YSwgZGF5dGltZT09JzAnJiBET1k+MzA1IHwgRE9ZIDwxMjEgKSRjbzJfZmx1eCwgXG4gICAgIFxuICAgICAgICAgICAgICAgeWxhYj1leHByZXNzaW9uKEZDficoJ35tdX5tb2x+bV57LTJ9fnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIHhsYWI9IGV4cHJlc3Npb24odVsnKiddficoJ35tfnNeey0xfX4nKScpLFxuICAgICB5bGltID0gYygwLDEwKSxcbiAgICAgICAgICAgICAgIG1haW49J05pZ2h0dGltZSBXZXQgR3Jvd2luZyBTZWFzb24gKE5vdiAyLSBBcHJpbCAzMCkgJyxcbiAgICAgIGNleD0gMC41KVxuXG5cblxuXG4jIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyNcbiNOaWdodHRpbWUgcGxvdCB0ZW1wZXJhdHV0cmUgMC01IENcbnNjYXR0ZXIuc21vb3RoKHN1YnNldChjZGF0YSwgZGF5dGltZT09JzAnJiBhaXJfdGVtcGVyYXR1cmVfYWRqPj0wICYgYWlyX3RlbXBlcmF0dXJlX2FkaiA8PTUpJHUuLFxuICAgICBzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgYWlyX3RlbXBlcmF0dXJlX2Fkaj4wICYgYWlyX3RlbXBlcmF0dXJlX2FkaiA8PTUgKSRjbzJfZmx1eCwgXG4gICAgIFxuICAgICAgICAgICAgICAgeWxhYj1leHByZXNzaW9uKEZDficoJ35tdX5tb2x+bV57LTJ9fnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIHhsYWI9IGV4cHJlc3Npb24odVsnKiddficoJ35tfnNeey0xfX4nKScpLFxuICAgICB5bGltID0gYygwLDEwKSxcbiAgICAgICAgICAgICAgIG1haW49J05pZ2h0dGltZSAwLTUgQyAnLFxuICAgICAgY2V4PSAwLjUpXG5cbiNOaWdodHRpbWUgcGxvdCB0ZW1wZXJhdHV0cmUgNS0xMCBDXG5zY2F0dGVyLnNtb290aChzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgYWlyX3RlbXBlcmF0dXJlX2Fkaj49NSAmIGFpcl90ZW1wZXJhdHVyZV9hZGogPD0xMCkkdS4sXG4gICAgIHN1YnNldChjZGF0YSwgZGF5dGltZT09JzAnJiBhaXJfdGVtcGVyYXR1cmVfYWRqPj01ICYgYWlyX3RlbXBlcmF0dXJlX2FkaiA8PTEwICkkY28yX2ZsdXgsIFxuICAgICBcbiAgICAgICAgICAgICAgIHlsYWI9ZXhwcmVzc2lvbihGQ34nKCd+bXV+bW9sfm1eey0yfX5zXnstMX1+JyknKSwgXG4gICAgICAgICAgICAgICB4bGFiPSBleHByZXNzaW9uKHVbJyonXX4nKCd+bX5zXnstMX1+JyknKSxcbiAgICAgeWxpbSA9IGMoMCwxMCksXG4gICAgICAgICAgICAgICBtYWluPSdOaWdodHRpbWUgNS0xMCBDICcsXG4gICAgICBjZXg9IDAuNSlcblxuI05pZ2h0dGltZSBwbG90IHRlbXBlcmF0dXRyZSAxMC0xNSBDXG5zY2F0dGVyLnNtb290aChzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgYWlyX3RlbXBlcmF0dXJlX2Fkaj49MTAgJiBhaXJfdGVtcGVyYXR1cmVfYWRqIDw9MTUpJHUuLFxuICAgICBzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgYWlyX3RlbXBlcmF0dXJlX2Fkaj49MTAgJiBhaXJfdGVtcGVyYXR1cmVfYWRqIDw9MTUgKSRjbzJfZmx1eCwgXG4gICAgIFxuICAgICAgICAgICAgICAgeWxhYj1leHByZXNzaW9uKEZDficoJ35tdX5tb2x+bV57LTJ9fnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIHhsYWI9IGV4cHJlc3Npb24odVsnKiddficoJ35tfnNeey0xfX4nKScpLFxuICAgICB5bGltID0gYygwLDEwKSxcbiAgICAgICAgICAgICAgIG1haW49J05pZ2h0dGltZSAxMC0xNSBDICcsXG4gICAgICBjZXg9IDAuNSlcblxuI05pZ2h0dGltZSBwbG90IHRlbXBlcmF0dXRyZSAxNS0yMCBDXG5zY2F0dGVyLnNtb290aChzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgYWlyX3RlbXBlcmF0dXJlX2Fkaj49MTUgJiBhaXJfdGVtcGVyYXR1cmVfYWRqIDw9MjApJHUuLFxuICAgICBzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgYWlyX3RlbXBlcmF0dXJlX2Fkaj49MTUgJiBhaXJfdGVtcGVyYXR1cmVfYWRqIDw9MjAgKSRjbzJfZmx1eCwgXG4gICAgIFxuICAgICAgICAgICAgICAgeWxhYj1leHByZXNzaW9uKEZDficoJ35tdX5tb2x+bV57LTJ9fnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIHhsYWI9IGV4cHJlc3Npb24odVsnKiddficoJ35tfnNeey0xfX4nKScpLFxuICAgICB5bGltID0gYygwLDEwKSxcbiAgICAgICAgICAgICAgIG1haW49J05pZ2h0dGltZSAxNS0yMCBDICcsXG4gICAgICBjZXg9IDAuNSlcblxuXG4jTmlnaHR0aW1lIHBsb3QgdGVtcGVyYXR1dHJlIDIwLTI1IENcbnNjYXR0ZXIuc21vb3RoKHN1YnNldChjZGF0YSwgZGF5dGltZT09JzAnJiBhaXJfdGVtcGVyYXR1cmVfYWRqPj0yMCAmIGFpcl90ZW1wZXJhdHVyZV9hZGogPD0yNSkkdS4sXG4gICAgIHN1YnNldChjZGF0YSwgZGF5dGltZT09JzAnJiBhaXJfdGVtcGVyYXR1cmVfYWRqPj0yMCAmIGFpcl90ZW1wZXJhdHVyZV9hZGogPD0yNSApJGNvMl9mbHV4LCBcbiAgICAgXG4gICAgICAgICAgICAgICB5bGFiPWV4cHJlc3Npb24oRkN+Jygnfm11fm1vbH5tXnstMn1+c157LTF9ficpJyksIFxuICAgICAgICAgICAgICAgeGxhYj0gZXhwcmVzc2lvbih1WycqJ11+Jygnfm1+c157LTF9ficpJyksXG4gICAgIHlsaW0gPSBjKDAsMTApLFxuICAgICAgICAgICAgICAgbWFpbj0nTmlnaHR0aW1lIDIwLTI1IEMgJyxcbiAgICAgIGNleD0gMC41KVxuXG4jTmlnaHR0aW1lIHBsb3QgdGVtcGVyYXR1dHJlIDI1LTMwIENcbnNjYXR0ZXIuc21vb3RoKHN1YnNldChjZGF0YSwgZGF5dGltZT09JzAnJiBhaXJfdGVtcGVyYXR1cmVfYWRqPj0yNSAmIGFpcl90ZW1wZXJhdHVyZV9hZGogPD0zMCkkdS4sXG4gICAgIHN1YnNldChjZGF0YSwgZGF5dGltZT09JzAnJiBhaXJfdGVtcGVyYXR1cmVfYWRqPj0yNSAmIGFpcl90ZW1wZXJhdHVyZV9hZGogPD0zMCApJGNvMl9mbHV4LCBcbiAgICAgXG4gICAgICAgICAgICAgICB5bGFiPWV4cHJlc3Npb24oRkN+Jygnfm11fm1vbH5tXnstMn1+c157LTF9ficpJyksIFxuICAgICAgICAgICAgICAgeGxhYj0gZXhwcmVzc2lvbih1WycqJ11+Jygnfm1+c157LTF9ficpJyksXG4gICAgIHlsaW0gPSBjKDAsMTApLFxuICAgICAgICAgICAgICAgbWFpbj0nTmlnaHR0aW1lIDI1LTMwIEMgJyxcbiAgICAgIGNleD0gMC41KVxuXG4jTmlnaHR0aW1lIHBsb3QgdGVtcGVyYXR1dHJlIGdyZWF0ZXIgdGhhbiAzMCBDXG5zY2F0dGVyLnNtb290aChzdWJzZXQoY2RhdGEsIGRheXRpbWU9PScwJyYgYWlyX3RlbXBlcmF0dXJlX2Fkaj49MzAgKSR1LixcbiAgICAgc3Vic2V0KGNkYXRhLCBkYXl0aW1lPT0nMCcmIGFpcl90ZW1wZXJhdHVyZV9hZGo+PTMwICApJGNvMl9mbHV4LCBcbiAgICAgXG4gICAgICAgICAgICAgICB5bGFiPWV4cHJlc3Npb24oRkN+Jygnfm11fm1vbH5tXnstMn1+c157LTF9ficpJyksIFxuICAgICAgICAgICAgICAgeGxhYj0gZXhwcmVzc2lvbih1WycqJ11+Jygnfm1+c157LTF9ficpJyksXG4gICAgIHlsaW0gPSBjKDAsMTApLFxuICAgICAgICAgICAgICAgbWFpbj0nTmlnaHR0aW1lIEdyZWF0ZXIgdGhhbiAzMCBDICcsXG4gICAgICBjZXg9IDAuNSlcblxuXG5cbiMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjI1xuI3UqIGJ5IHRlbXBlcmF0dXJlXG5cbnNjYXR0ZXIuc21vb3RoKGNkYXRhJHUuLCBcbiAgICAgIGNkYXRhJGFpcl90ZW1wZXJhdHVyZV9hZGosXG4gICAgICAgICAgICAgICB4bGFiPWV4cHJlc3Npb24odVsnKiddficoJ35tfnNeey0xfX4nKScpLCBcbiAgICAgICAgICAgICAgIHlsYWI9IGV4cHJlc3Npb24oQWlyfnRlbXBlcmF0dXJlficoJ35kZWdyZWV+Q34nKScpLCBcbiAgICAgICAgICAgICAgIG1haW49JyAnKVxuXG4jbmlnaHR0aW1lIGxpbmVhciByZWxhdGlvbnNoaXAgYmV0d2VlbiBjbzIgZmx1eCBhbmQgVSZcbmxtX2Fpcl90ZW1wX3U8LSBsbSh1LiB+XG4gICAgICAgICAgICAgICAgICAgICBhaXJfdGVtcGVyYXR1cmVfYWRqXG4sIGRhdGEgPSBjZGF0YSlcblxuc3VtbWFyeShsbV9haXJfdGVtcF91KVxuXG5cblxuYGBgXG5gYGAifQ== -->
```r
```r
######################################################################
#Daytime and nighttime plot
scatter.smooth( cdata$u., cdata$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
main='Daytime and Nighttime')
#day and night linear relationship
lm_all_time<- lm(u. ~
co2_flux
,
data = cdata)
summary(lm_all_time)
#############################################
#Daytime plot
scatter.smooth(subset(cdata, daytime=='1')$u.,
subset(cdata, daytime=='1')$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
main='Daytime ')
#####################################################################
#daytime linear relationship between co2 flux and U&
lm_day_time<- lm(subset(cdata, daytime=='1')$u. ~
subset(cdata, daytime=='1')$co2_flux
,
data = cdata)
summary(lm_day_time)
########################################################
#Nighttime plot feb_march
smoothScatter(subset(cdata, daytime=='0'& DOY>31 & DOY <120)$u.,
subset(cdata, daytime=='0'& DOY>31 & DOY <120 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
main='Nighttime ',
cex= 0.5)
#nighttime linear relationship between co2 flux and U&
lm_night_time_feb_mar<- lm(subset(cdata, daytime=='0'& DOY>31 & DOY <120)$u.~
subset(cdata, daytime=='0'& DOY>31 & DOY <120 )$co2_flux
, data = cdata)
summary(lm_night_time_feb_mar)
#######################################
#Nighttime plot Dry Season. May 1(121) -Nov 1(305)
scatter.smooth(subset(cdata, daytime=='0'& DOY>120 & DOY <306)$u.,
subset(cdata, daytime=='0'& DOY>120 & DOY <306 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime Dry Season (May 1- Nov. 1) ',
cex= 0.5)
##############################
#Nighttime plot Wet Season, Nov 2- April 30
scatter.smooth(subset(cdata, daytime=='0'& DOY>305 | DOY <121)$u.,
subset(cdata, daytime=='0'& DOY>305 | DOY <121 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime Wet Growing Season (Nov 2- April 30) ',
cex= 0.5)
#######################################
#Nighttime plot temperatutre 0-5 C
scatter.smooth(subset(cdata, daytime=='0'& air_temperature_adj>=0 & air_temperature_adj <=5)$u.,
subset(cdata, daytime=='0'& air_temperature_adj>0 & air_temperature_adj <=5 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime 0-5 C ',
cex= 0.5)
#Nighttime plot temperatutre 5-10 C
scatter.smooth(subset(cdata, daytime=='0'& air_temperature_adj>=5 & air_temperature_adj <=10)$u.,
subset(cdata, daytime=='0'& air_temperature_adj>=5 & air_temperature_adj <=10 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime 5-10 C ',
cex= 0.5)
#Nighttime plot temperatutre 10-15 C
scatter.smooth(subset(cdata, daytime=='0'& air_temperature_adj>=10 & air_temperature_adj <=15)$u.,
subset(cdata, daytime=='0'& air_temperature_adj>=10 & air_temperature_adj <=15 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime 10-15 C ',
cex= 0.5)
#Nighttime plot temperatutre 15-20 C
scatter.smooth(subset(cdata, daytime=='0'& air_temperature_adj>=15 & air_temperature_adj <=20)$u.,
subset(cdata, daytime=='0'& air_temperature_adj>=15 & air_temperature_adj <=20 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime 15-20 C ',
cex= 0.5)
#Nighttime plot temperatutre 20-25 C
scatter.smooth(subset(cdata, daytime=='0'& air_temperature_adj>=20 & air_temperature_adj <=25)$u.,
subset(cdata, daytime=='0'& air_temperature_adj>=20 & air_temperature_adj <=25 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime 20-25 C ',
cex= 0.5)
#Nighttime plot temperatutre 25-30 C
scatter.smooth(subset(cdata, daytime=='0'& air_temperature_adj>=25 & air_temperature_adj <=30)$u.,
subset(cdata, daytime=='0'& air_temperature_adj>=25 & air_temperature_adj <=30 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime 25-30 C ',
cex= 0.5)
#Nighttime plot temperatutre greater than 30 C
scatter.smooth(subset(cdata, daytime=='0'& air_temperature_adj>=30 )$u.,
subset(cdata, daytime=='0'& air_temperature_adj>=30 )$co2_flux,
ylab=expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
xlab= expression(u['*']~'('~m~s^{-1}~')'),
ylim = c(0,10),
main='Nighttime Greater than 30 C ',
cex= 0.5)
################################################
#u* by temperature
scatter.smooth(cdata$u.,
cdata$air_temperature_adj,
xlab=expression(u['*']~'('~m~s^{-1}~')'),
ylab= expression(Air~temperature~'('~degree~C~')'),
main=' ')
#nighttime linear relationship between co2 flux and U&
lm_air_temp_u<- lm(u. ~
air_temperature_adj
, data = cdata)
summary(lm_air_temp_u)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#diurnal pattern of soil heat flux plate 1. Also did 2 to see if it is correct. SHF1 is wack for June.
#Looks like it began messing up in May
#plot(tapply(combo_master_ed_met$LE ,round(combo_master_ed_met$DOY),function(x) mean(x,na.rm=T))
#subset(cdata, daytime=='0'& DOY>120 & DOY <306)$u
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuXG4jYWxsIHRpbWUgcGF0dGVybiBvZiBTSEYxIHBsYXRlXG5wbG90KHRhcHBseShjZGF0YSRDb3JyZWN0X3NoZl8xLCBjZGF0YSR0aW1lLGZ1bmN0aW9uKHgpIG1lYW4oeCxuYS5ybT1UKSksIG1haW49IFxcU0hGMSBBbGx0aW1lXFwpXG5cbiNkaXVyaW5hbCBwYXR0ZXJuIFNIRjEgZm9yIEFwcmlsLiAjIGxvb2tzIGdvb2RcbnBsb3QodGFwcGx5KHN1YnNldChjZGF0YSwgRE9ZID49IDkyICYgRE9ZPD0gMTIxKSRDb3JyZWN0X3NoZl8xLCBzdWJzZXQoY2RhdGEsIERPWSA+PSA5MiAmIERPWTw9IDEyMSkgJHRpbWUsZnVuY3Rpb24oeCkgbWVhbih4LG5hLnJtPVQpKSwgbWFpbj0gXFxTSEYxIEFwcmlsXFwpXG5cblxuI2RpdXJpbmFsIHBhdHRlcm4gU0hGMSBmb3IgTWF5LiAgbG9va3Mgd2Fjay4gYnV0IG5vdCBhcyBiYWQgYXMgSnVuZVxucGxvdCh0YXBwbHkoc3Vic2V0KGNkYXRhLCBET1kgPj0gMTIyICYgRE9ZPD0gMTUyKSRDb3JyZWN0X3NoZl8xLCBzdWJzZXQoY2RhdGEsIERPWSA+PSAxMjIgJiBET1k8PSAxNTIpICR0aW1lLGZ1bmN0aW9uKHgpIG1lYW4oeCxuYS5ybT1UKSksIG1haW49IFxcU0hGMSBNYXlcXClcblxucGxvdCh0YXBwbHkoc3Vic2V0KGNkYXRhLCBUSU1FU1RBTVA+PSBcXDIwMjAtMDctMjIgMDowMCBHTVQtOFxcICkkQ29ycmVjdF9zaGZfMSwgc3Vic2V0KGNkYXRhLCBUSU1FU1RBTVA+PSBcXDIwMjAtMDctMjIgMDowMCBHTVQtOFxcICkkdGltZSxmdW5jdGlvbih4KSBtZWFuKHgsbmEucm09VCkpLCBtYWluPSBcXFNIRjE+IEp1bHkgMjFcXCwgeGxhYiA9IFxcRXZlcnkgMzAgbWludXRlc1xcLCB5bGFiID0gXFxTb2lsIEhlYXQgRmx1eFxcKVxuXG5wbG90X0p1bHlfMjAyMDwtcGxvdCh0YXBwbHkoc3Vic2V0KGNkYXRhLCBUSU1FU1RBTVA+PSBcXDIwMjAtMDctMjIgMDowMCBHTVQtOFxcICkkQ29ycmVjdF9zaGZfMiwgc3Vic2V0KGNkYXRhLCBUSU1FU1RBTVA+PSBcXDIwMjAtMDctMjIgMDowMCBHTVQtOFxcICkkdGltZSxmdW5jdGlvbih4KSBtZWFuKHgsbmEucm09VCkpLCBtYWluPSBcXFNIRjI+IEp1bHkgMjFcXCwgeGxhYiA9IFxcRXZlcnkgMzAgbWludXRlc1xcLCB5bGFiID0gXFxTb2lsIEhlYXQgRmx1eFxcKVxuXG5cblxuXG4jY2RhdGEkVElNRVNUQU1QJG1vbjw9M3xjZGF0YSRUSU1FU1RBTVAkbW9uPj0xMVxuIyBmb3JtYXQ9XFwlbS8lZC8lWSAlSDolTVxcLCAgdHogPSBcXEV0Yy9HTVQtOFxcKVxuXG4jICNkaXVyaW5hbCBwYXR0ZXJuIFNIRjEgZm9yIEp1bmUuICMgbG9va3Mgd2FjayBzdGlsbCBmaWx0ZXJlZCBvdXRcbiMgcGxvdCh0YXBwbHkoc3Vic2V0KGNkYXRhLCBET1kgPj0gMTUzICYgRE9ZPD0gMTc4KSRDb3JyZWN0X3NoZl8xLCBzdWJzZXQoY2RhdGEsIERPWSA+PSAxNTMgJiBET1k8PSAxNzgpICR0aW1lLGZ1bmN0aW9uKHgpIG1lYW4oeCxuYS5ybT1UKSksIG1haW49IFxcU0hGMSBKdW5lXFwpXG5cbiMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjI1xuI1NIRjIgZm9yIGNvbXBhcmlzb25cblxuI2FsbCB0aW1lIHBhdHRlcm4gb2YgU0hGMSBwbGF0ZVxucGxvdCh0YXBwbHkoY2RhdGEkQ29ycmVjdF9zaGZfMiwgY2RhdGEkdGltZSxmdW5jdGlvbih4KSBtZWFuKHgsbmEucm09VCkpLCBtYWluPSBcXFNIRjIgQWxsIFRpbWVcXClcblxuI2RpdXJpbmFsIHBhdHRlcm4gU0hGMiBmb3IgQXByaWwuICMgbG9va3MgZ29vZFxucGxvdCh0YXBwbHkoc3Vic2V0KGNkYXRhLCBET1kgPj0gOTIgJiBET1k8PSAxMjEpJENvcnJlY3Rfc2hmXzIsIHN1YnNldChjZGF0YSwgRE9ZID49IDkyICYgRE9ZPD0gMTIxKSAkdGltZSxmdW5jdGlvbih4KSBtZWFuKHgsbmEucm09VCkpLCBtYWluPSBcXFNIRjIgQXByaWxcXClcblxuXG4jZGl1cmluYWwgcGF0dGVybiBTSEYyIGZvciBNYXkuIFxucGxvdCh0YXBwbHkoc3Vic2V0KGNkYXRhLCBET1kgPj0gMTIyICYgRE9ZPD0gMTUyKSRDb3JyZWN0X3NoZl8yLCBzdWJzZXQoY2RhdGEsIERPWSA+PSAxMjIgJiBET1k8PSAxNTIpICR0aW1lLGZ1bmN0aW9uKHgpIG1lYW4oeCxuYS5ybT1UKSksIG1haW49IFxcU0hGMiBNYXlcXClcblxuXG5cbiNkaXVyaW5hbCBwYXR0ZXJuIFNIRjIgZm9yIEp1bmUuICMgbG9va3MgZ29vZC4gU28gU0hGMSBpcyBlbWl0dGluZyBub2lzZS4gXG5wbG90KHRhcHBseShzdWJzZXQoY2RhdGEsIERPWSA+PSAxNTMgJiBET1k8PSAxNzgpJENvcnJlY3Rfc2hmXzIsIHN1YnNldChjZGF0YSwgRE9ZID49IDE1MyAmIERPWTw9IDE3OCkgJHRpbWUsZnVuY3Rpb24oeCkgbWVhbih4LG5hLnJtPVQpKSwgbWFpbj0gXFxTSEYgMiBKdW5lXFwpXG5cblxuYGBgXG5gYGAifQ== -->
```r
```r
#all time pattern of SHF1 plate
plot(tapply(cdata$Correct_shf_1, cdata$time,function(x) mean(x,na.rm=T)), main= \SHF1 Alltime\)
#diurinal pattern SHF1 for April. # looks good
plot(tapply(subset(cdata, DOY >= 92 & DOY<= 121)$Correct_shf_1, subset(cdata, DOY >= 92 & DOY<= 121) $time,function(x) mean(x,na.rm=T)), main= \SHF1 April\)
#diurinal pattern SHF1 for May. looks wack. but not as bad as June
plot(tapply(subset(cdata, DOY >= 122 & DOY<= 152)$Correct_shf_1, subset(cdata, DOY >= 122 & DOY<= 152) $time,function(x) mean(x,na.rm=T)), main= \SHF1 May\)
plot(tapply(subset(cdata, TIMESTAMP>= \2020-07-22 0:00 GMT-8\ )$Correct_shf_1, subset(cdata, TIMESTAMP>= \2020-07-22 0:00 GMT-8\ )$time,function(x) mean(x,na.rm=T)), main= \SHF1> July 21\, xlab = \Every 30 minutes\, ylab = \Soil Heat Flux\)
plot_July_2020<-plot(tapply(subset(cdata, TIMESTAMP>= \2020-07-22 0:00 GMT-8\ )$Correct_shf_2, subset(cdata, TIMESTAMP>= \2020-07-22 0:00 GMT-8\ )$time,function(x) mean(x,na.rm=T)), main= \SHF2> July 21\, xlab = \Every 30 minutes\, ylab = \Soil Heat Flux\)
#cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11
# format=\%m/%d/%Y %H:%M\, tz = \Etc/GMT-8\)
# #diurinal pattern SHF1 for June. # looks wack still filtered out
# plot(tapply(subset(cdata, DOY >= 153 & DOY<= 178)$Correct_shf_1, subset(cdata, DOY >= 153 & DOY<= 178) $time,function(x) mean(x,na.rm=T)), main= \SHF1 June\)
####################################
#SHF2 for comparison
#all time pattern of SHF1 plate
plot(tapply(cdata$Correct_shf_2, cdata$time,function(x) mean(x,na.rm=T)), main= \SHF2 All Time\)
#diurinal pattern SHF2 for April. # looks good
plot(tapply(subset(cdata, DOY >= 92 & DOY<= 121)$Correct_shf_2, subset(cdata, DOY >= 92 & DOY<= 121) $time,function(x) mean(x,na.rm=T)), main= \SHF2 April\)
#diurinal pattern SHF2 for May.
plot(tapply(subset(cdata, DOY >= 122 & DOY<= 152)$Correct_shf_2, subset(cdata, DOY >= 122 & DOY<= 152) $time,function(x) mean(x,na.rm=T)), main= \SHF2 May\)
#diurinal pattern SHF2 for June. # looks good. So SHF1 is emitting noise.
plot(tapply(subset(cdata, DOY >= 153 & DOY<= 178)$Correct_shf_2, subset(cdata, DOY >= 153 & DOY<= 178) $time,function(x) mean(x,na.rm=T)), main= \SHF 2 June\)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
# Sensible Heat, Latent Heat U* and Co2 fluxes by Net radiation, seperated by wind directions
#Treatment 150-230, Control 230-310
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuV0QuU0VfTEVfY28yX2J5X05SPC1yZXAoMixucm93KGNkYXRhKSlcbldELlNFX0xFX2NvMl9ieV9OUlt3aGljaCghaXMubmEoY2RhdGEkd2luZF9kaXIpJlxuICAgICAgICAgICAgICAgKGNkYXRhJHdpbmRfZGlyPjE1MCAmIGNkYXRhJHdpbmRfZGlyPD0yMzApJiAjdHJlYXRtZW50XG4gICAgICAgICAgICAgICAgIWlzLm5hKGNkYXRhJHdpbmRfZGlyKSxcbiAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2Rpcj4yMzAgJiBjZGF0YSR3aW5kX2Rpcjw9MzEwKV08LTEgI2NvbnRyb2xcbldkLmxlZ2VuZF8xPC1jKFxcMTUwLTIzMCB3aW5kIFxcLFxcMjMwLTMxMCB3aW5kXFwpXG5XZC5jb2xfMTwtYyhyZ2IoMCwwLDEsMC41LG1heENvbG9yVmFsdWU9MSkscmdiKDEsMCwwLDAuNSxtYXhDb2xvclZhbHVlPTEpKVxuXG50YXJnZXQucGxvdC52YXJfMTwtYyhcXGNvMl9mbHV4XFwsXFxMRVxcLFxcSFxcLFxcdS5cXClcbnRhcmdldC5wbG90LnZhcl8xLnRpdGxlPC1jKGV4cHJlc3Npb24oRkN+Jygnfm11fm1vbH5tXnstMn1+c157LTF9ficpJyksXG4gICAgICAgICAgICAgICAgICAgICAgICAgZXhwcmVzc2lvbihMRX4nKCd+V35tXnstMn1+JyknKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICBleHByZXNzaW9uKEh+Jygnfld+bV57LTJ9ficpJyksXG4gICAgICAgICAgICAgICAgICAgICAgICAgZXhwcmVzc2lvbih1WycqJ11+Jygnfm1+c157LTF9ficpJykpIFxuXG5mb3IoazEgaW4gMTpsZW5ndGgodGFyZ2V0LnBsb3QudmFyXzEpKXtcbiAgXG4gICMjIGxvY2F0ZSB0aGUgc3RhcnQgb2YgZWFjaCBtb250aFxuICBtb250aC5sb2MxPC13aGljaChjZGF0YSRUSU1FU1RBTVAkbWRheT09MSZcbiAgICAgICAgICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCRob3VyPT0wJlxuICAgICAgICAgICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJG1pbj09MClcbiAgbW9udGgudGlja3MxIDwtIHNlcShjZGF0YSRUSU1FU1RBTVBbbW9udGgubG9jMVsxXV0sXG4gICAgICAgICAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVBbbW9udGgubG9jMVtsZW5ndGgobW9udGgubG9jMSldXSxieT1cXG1vbnRoc1xcKVxuICBcbiAgcG5nKHBhc3RlMChwbG90LnBhdGgsXG4gICAgICAgICAgICAgXFxDb25jb3JkX0hfTEVfY28yXzFfXFwsXG4gICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJHllYXJbMV0rMTkwMCxcXF9cXCxcbiAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVAkeWRheVsxXSsxLFxcX1xcLFxuICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCR5ZWFyW25yb3coY2RhdGEpXSsxOTAwLFxcX1xcLFxuICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCR5ZGF5W25yb3coY2RhdGEpXSsxLFxcX1xcLFxuICAgICAgICAgICAgIHRhcmdldC5wbG90LnZhcl8xW2sxXSxcXF9jb2xvcl9cXCxcbiAgICAgICAgICAgICBTeXMuRGF0ZSgpLFxcLnBuZ1xcKSxcbiAgICAgIHdpZHRoPTgsXG4gICAgICBoZWlnaHQ9NCxcbiAgICAgIHVuaXRzPVxcaW5cXCxcbiAgICAgIHJlcz0zMDAsXG4gICAgICBwb2ludHNpemUgPSAxMSxcbiAgICAgIGJnID0gXFx3aGl0ZVxcKVxuICBcbiAgcGFyKG9tYT1jKDQsNC41LDAuNSwwLjUpLG1hcj1jKDAsMC41LDAsMC41KSxmaWc9YygwLDAuNywwLDEpKVxuICBwbG90KGNkYXRhJENvcnJlY3RfTlIsXG4gICAgICAgY2RhdGFbLHRhcmdldC5wbG90LnZhcl8xW2sxXV0sXG4gICAgICAgeGxhYj1cXFxcLFxuICAgICAgIHlsYWI9XFxcXCxcbiAgICAgICBjZXg9MC41LGNvbD1XZC5jb2xfMVtXRC5TRV9MRV9jbzJfYnlfTlJdLFxuICAgICAgIGxhcz0xLHBjaD0xNixcbiAgICAgICB4YXhzPVxcaVxcLHlheHM9XFxpXFxcbiAgKVxuICBtdGV4dChzaWRlPTIsdGFyZ2V0LnBsb3QudmFyXzEudGl0bGVbazFdLGxpbmU9MylcbiAgbXRleHQoc2lkZT0xLFxcTmV0IFJhZGlhdGlvblxcLGxpbmU9Mi44KVxuICBhYmxpbmUoaD0wLGNvbD1cXGRhcmtncmV5XFwpXG4gIGF4aXMoMSwgYXQgPSBtb250aC50aWNrczEsIGxhYmVscyA9IEZBTFNFLCB0Y2wgPSAtMC4zKVxuICBcbiAgcGFyKGZpZz1jKDAuNywxLDAsMSksbmV3PVQpXG4gIGhpc3QwQTwtaGlzdChjZGF0YVssdGFyZ2V0LnBsb3QudmFyXzFbazFdXSxcbiAgICAgICAgICAgICAgcGxvdD1GLG5jbGFzcz01MClcbiAgaGlzdDFBPC1oaXN0KGNkYXRhWyx0YXJnZXQucGxvdC52YXJfMVtrMV1dW1dELlNFX0xFX2NvMl9ieV9OUj09MV0sXG4gICAgICAgICAgICAgIHBsb3Q9RixicmVha3M9aGlzdDBBJGJyZWFrcylcbiAgaGlzdDJBPC1oaXN0KGNkYXRhWyx0YXJnZXQucGxvdC52YXJfMVtrMV1dW1dELlNFX0xFX2NvMl9ieV9OUj09Ml0sXG4gICAgICAgICAgICAgIHBsb3Q9RixicmVha3M9aGlzdDBBJGJyZWFrcykgXG4gIFxuICBiYXJwbG90KGhpc3QxQSRjb3VudHMsXG4gICAgICAgICAgYXhlcz1GLFxuICAgICAgICAgIGhvcml6PVQsXG4gICAgICAgICAgeWxpbT1jKDAsbGVuZ3RoKGhpc3QxQSRicmVha3MpKzEpLFxuICAgICAgICAgIHhsaW09YygtNSxtYXgoYyhoaXN0MUEkY291bnRzLGhpc3QyQSRjb3VudHMpKSksXG4gICAgICAgICAgc3BhY2U9MCxjb2w9V2QuY29sXzFbMV0sYm9yZGVyPU5BKSAjIGJhcnBsb3RcbiAgYmFycGxvdChoaXN0MkEkY291bnRzLFxuICAgICAgICAgIGF4ZXM9RixcbiAgICAgICAgICBhZGQ9VCxcbiAgICAgICAgICBob3Jpej1ULFxuICAgICAgICAgIHlsaW09YygwLGxlbmd0aChoaXN0MUEkYnJlYWtzKSsxKSxcbiAgICAgICAgICB4bGltPWMoLTUsbWF4KGMoaGlzdDFBJGNvdW50cyxoaXN0MkEkY291bnRzKSkpLFxuICAgICAgICAgIHNwYWNlPTAsY29sPVdkLmNvbF8xWzJdLGJvcmRlcj1OQSkgIyBiYXJwbG90XG4gIGxlZ2VuZCgwLFxuICAgICAgICAgbGVuZ3RoKGhpc3QxQSRicmVha3MpKzEsXG4gICAgICAgICBmaWxsPVdkLmNvbF8xLGJvcmRlcj1OQSxcbiAgICAgICAgIGxlZ2VuZD1XZC5sZWdlbmRfMSxidHk9XFxuXFwsXG4gICAgICAgICBjZXg9MC45KVxuICBkZXYub2ZmKClcbn1cbiAgXG5gYGBcbmBgYCJ9 -->
```r
```r
WD.SE_LE_co2_by_NR<-rep(2,nrow(cdata))
WD.SE_LE_co2_by_NR[which(!is.na(cdata$wind_dir)&
(cdata$wind_dir>150 & cdata$wind_dir<=230)& #treatment
!is.na(cdata$wind_dir),
cdata$wind_dir>230 & cdata$wind_dir<=310)]<-1 #control
Wd.legend_1<-c(\150-230 wind \,\230-310 wind\)
Wd.col_1<-c(rgb(0,0,1,0.5,maxColorValue=1),rgb(1,0,0,0.5,maxColorValue=1))
target.plot.var_1<-c(\co2_flux\,\LE\,\H\,\u.\)
target.plot.var_1.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(LE~'('~W~m^{-2}~')'),
expression(H~'('~W~m^{-2}~')'),
expression(u['*']~'('~m~s^{-1}~')'))
for(k1 in 1:length(target.plot.var_1)){
## locate the start of each month
month.loc1<-which(cdata$TIMESTAMP$mday==1&
cdata$TIMESTAMP$hour==0&
cdata$TIMESTAMP$min==0)
month.ticks1 <- seq(cdata$TIMESTAMP[month.loc1[1]],
cdata$TIMESTAMP[month.loc1[length(month.loc1)]],by=\months\)
png(paste0(plot.path,
\Concord_H_LE_co2_1_\,
cdata$TIMESTAMP$year[1]+1900,\_\,
cdata$TIMESTAMP$yday[1]+1,\_\,
cdata$TIMESTAMP$year[nrow(cdata)]+1900,\_\,
cdata$TIMESTAMP$yday[nrow(cdata)]+1,\_\,
target.plot.var_1[k1],\_color_\,
Sys.Date(),\.png\),
width=8,
height=4,
units=\in\,
res=300,
pointsize = 11,
bg = \white\)
par(oma=c(4,4.5,0.5,0.5),mar=c(0,0.5,0,0.5),fig=c(0,0.7,0,1))
plot(cdata$Correct_NR,
cdata[,target.plot.var_1[k1]],
xlab=\\,
ylab=\\,
cex=0.5,col=Wd.col_1[WD.SE_LE_co2_by_NR],
las=1,pch=16,
xaxs=\i\,yaxs=\i\
)
mtext(side=2,target.plot.var_1.title[k1],line=3)
mtext(side=1,\Net Radiation\,line=2.8)
abline(h=0,col=\darkgrey\)
axis(1, at = month.ticks1, labels = FALSE, tcl = -0.3)
par(fig=c(0.7,1,0,1),new=T)
hist0A<-hist(cdata[,target.plot.var_1[k1]],
plot=F,nclass=50)
hist1A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==1],
plot=F,breaks=hist0A$breaks)
hist2A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==2],
plot=F,breaks=hist0A$breaks)
barplot(hist1A$counts,
axes=F,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[1],border=NA) # barplot
barplot(hist2A$counts,
axes=F,
add=T,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[2],border=NA) # barplot
legend(0,
length(hist1A$breaks)+1,
fill=Wd.col_1,border=NA,
legend=Wd.legend_1,bty=\n\,
cex=0.9)
dev.off()
}
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
# Sensible Heat, Latent Heat U* and Co2 fluxes by Net radiation, seperated by wind directions
#Treatment 150-230, Control 230-280
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuV0QuU0VfTEVfY28yX2J5X05SPC1yZXAoMixucm93KGNkYXRhKSlcbldELlNFX0xFX2NvMl9ieV9OUlt3aGljaCghaXMubmEoY2RhdGEkd2luZF9kaXIpJlxuICAgICAgICAgICAgICAgKGNkYXRhJHdpbmRfZGlyPjE1MCAmIGNkYXRhJHdpbmRfZGlyPD0yMzApJiAjdHJlYXRtZW50XG4gICAgICAgICAgICAgICAgIWlzLm5hKGNkYXRhJHdpbmRfZGlyKSxcbiAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2Rpcj4yMzAgJiBjZGF0YSR3aW5kX2Rpcjw9MjgwKV08LTEgI2NvbnRyb2xcbldkLmxlZ2VuZF8xPC1jKFwiMTUwLTIzMCB3aW5kIFwiLFwiMjMwLTI4MCB3aW5kXCIpXG5XZC5jb2xfMTwtYyhyZ2IoMCwwLDEsMC41LG1heENvbG9yVmFsdWU9MSkscmdiKDEsMCwwLDAuNSxtYXhDb2xvclZhbHVlPTEpKVxuXG50YXJnZXQucGxvdC52YXJfMTwtYyhcImNvMl9mbHV4XCIsXCJMRVwiLFwiSFwiLFwidS5cIilcbnRhcmdldC5wbG90LnZhcl8xLnRpdGxlPC1jKGV4cHJlc3Npb24oRkN+Jygnfm11fm1vbH5tXnstMn1+c157LTF9ficpJyksXG4gICAgICAgICAgICAgICAgICAgICAgICAgZXhwcmVzc2lvbihMRX4nKCd+V35tXnstMn1+JyknKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICBleHByZXNzaW9uKEh+Jygnfld+bV57LTJ9ficpJyksXG4gICAgICAgICAgICAgICAgICAgICAgICAgZXhwcmVzc2lvbih1WycqJ11+Jygnfm1+c157LTF9ficpJykpIFxuXG5mb3IoazEgaW4gMTpsZW5ndGgodGFyZ2V0LnBsb3QudmFyXzEpKXtcbiAgXG4gICMjIGxvY2F0ZSB0aGUgc3RhcnQgb2YgZWFjaCBtb250aFxuICBtb250aC5sb2MxPC13aGljaChjZGF0YSRUSU1FU1RBTVAkbWRheT09MSZcbiAgICAgICAgICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCRob3VyPT0wJlxuICAgICAgICAgICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJG1pbj09MClcbiAgbW9udGgudGlja3MxIDwtIHNlcShjZGF0YSRUSU1FU1RBTVBbbW9udGgubG9jMVsxXV0sXG4gICAgICAgICAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVBbbW9udGgubG9jMVtsZW5ndGgobW9udGgubG9jMSldXSxieT1cIm1vbnRoc1wiKVxuICBcbiAgcG5nKHBhc3RlMChwbG90LnBhdGgsXG4gICAgICAgICAgICAgXCJDb25jb3JkX0hfTEVfY28yXzJfXCIsXG4gICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJHllYXJbMV0rMTkwMCxcIl9cIixcbiAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVAkeWRheVsxXSsxLFwiX1wiLFxuICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCR5ZWFyW25yb3coY2RhdGEpXSsxOTAwLFwiX1wiLFxuICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCR5ZGF5W25yb3coY2RhdGEpXSsxLFwiX1wiLFxuICAgICAgICAgICAgIHRhcmdldC5wbG90LnZhcl8xW2sxXSxcIl9jb2xvcl9cIixcbiAgICAgICAgICAgICBTeXMuRGF0ZSgpLFwiLnBuZ1wiKSxcbiAgICAgIHdpZHRoPTgsXG4gICAgICBoZWlnaHQ9NCxcbiAgICAgIHVuaXRzPVwiaW5cIixcbiAgICAgIHJlcz0zMDAsXG4gICAgICBwb2ludHNpemUgPSAxMSxcbiAgICAgIGJnID0gXCJ3aGl0ZVwiKVxuICBcbiAgcGFyKG9tYT1jKDQsNC41LDAuNSwwLjUpLG1hcj1jKDAsMC41LDAsMC41KSxmaWc9YygwLDAuNywwLDEpKVxuICBwbG90KGNkYXRhJENvcnJlY3RfTlIsXG4gICAgICAgY2RhdGFbLHRhcmdldC5wbG90LnZhcl8xW2sxXV0sXG4gICAgICAgeGxhYj1cIlwiLFxuICAgICAgIHlsYWI9XCJcIixcbiAgICAgICBjZXg9MC41LGNvbD1XZC5jb2xfMVtXRC5TRV9MRV9jbzJfYnlfTlJdLFxuICAgICAgIGxhcz0xLHBjaD0xNixcbiAgICAgICB4YXhzPVwiaVwiLHlheHM9XCJpXCJcbiAgKVxuICBtdGV4dChzaWRlPTIsdGFyZ2V0LnBsb3QudmFyXzEudGl0bGVbazFdLGxpbmU9MylcbiAgbXRleHQoc2lkZT0xLFwiTmV0IFJhZGlhdGlvblwiLGxpbmU9Mi44KVxuICBhYmxpbmUoaD0wLGNvbD1cImRhcmtncmV5XCIpXG4gIGF4aXMoMSwgYXQgPSBtb250aC50aWNrczEsIGxhYmVscyA9IEZBTFNFLCB0Y2wgPSAtMC4zKVxuICBcbiAgcGFyKGZpZz1jKDAuNywxLDAsMSksbmV3PVQpXG4gIGhpc3QwQTwtaGlzdChjZGF0YVssdGFyZ2V0LnBsb3QudmFyXzFbazFdXSxcbiAgICAgICAgICAgICAgcGxvdD1GLG5jbGFzcz01MClcbiAgaGlzdDFBPC1oaXN0KGNkYXRhWyx0YXJnZXQucGxvdC52YXJfMVtrMV1dW1dELlNFX0xFX2NvMl9ieV9OUj09MV0sXG4gICAgICAgICAgICAgIHBsb3Q9RixicmVha3M9aGlzdDBBJGJyZWFrcylcbiAgaGlzdDJBPC1oaXN0KGNkYXRhWyx0YXJnZXQucGxvdC52YXJfMVtrMV1dW1dELlNFX0xFX2NvMl9ieV9OUj09Ml0sXG4gICAgICAgICAgICAgIHBsb3Q9RixicmVha3M9aGlzdDBBJGJyZWFrcykgXG4gIFxuICBiYXJwbG90KGhpc3QxQSRjb3VudHMsXG4gICAgICAgICAgYXhlcz1GLFxuICAgICAgICAgIGhvcml6PVQsXG4gICAgICAgICAgeWxpbT1jKDAsbGVuZ3RoKGhpc3QxQSRicmVha3MpKzEpLFxuICAgICAgICAgIHhsaW09YygtNSxtYXgoYyhoaXN0MUEkY291bnRzLGhpc3QyQSRjb3VudHMpKSksXG4gICAgICAgICAgc3BhY2U9MCxjb2w9V2QuY29sXzFbMV0sYm9yZGVyPU5BKSAjIGJhcnBsb3RcbiAgYmFycGxvdChoaXN0MkEkY291bnRzLFxuICAgICAgICAgIGF4ZXM9RixcbiAgICAgICAgICBhZGQ9VCxcbiAgICAgICAgICBob3Jpej1ULFxuICAgICAgICAgIHlsaW09YygwLGxlbmd0aChoaXN0MUEkYnJlYWtzKSsxKSxcbiAgICAgICAgICB4bGltPWMoLTUsbWF4KGMoaGlzdDFBJGNvdW50cyxoaXN0MkEkY291bnRzKSkpLFxuICAgICAgICAgIHNwYWNlPTAsY29sPVdkLmNvbF8xWzJdLGJvcmRlcj1OQSkgIyBiYXJwbG90XG4gIGxlZ2VuZCgwLFxuICAgICAgICAgbGVuZ3RoKGhpc3QxQSRicmVha3MpKzEsXG4gICAgICAgICBmaWxsPVdkLmNvbF8xLGJvcmRlcj1OQSxcbiAgICAgICAgIGxlZ2VuZD1XZC5sZWdlbmRfMSxidHk9XCJuXCIsXG4gICAgICAgICBjZXg9MC45KVxuICBkZXYub2ZmKClcbn1cbiAgXG5gYGAifQ== -->
```r
WD.SE_LE_co2_by_NR<-rep(2,nrow(cdata))
WD.SE_LE_co2_by_NR[which(!is.na(cdata$wind_dir)&
(cdata$wind_dir>150 & cdata$wind_dir<=230)& #treatment
!is.na(cdata$wind_dir),
cdata$wind_dir>230 & cdata$wind_dir<=280)]<-1 #control
Wd.legend_1<-c("150-230 wind ","230-280 wind")
Wd.col_1<-c(rgb(0,0,1,0.5,maxColorValue=1),rgb(1,0,0,0.5,maxColorValue=1))
target.plot.var_1<-c("co2_flux","LE","H","u.")
target.plot.var_1.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(LE~'('~W~m^{-2}~')'),
expression(H~'('~W~m^{-2}~')'),
expression(u['*']~'('~m~s^{-1}~')'))
for(k1 in 1:length(target.plot.var_1)){
## locate the start of each month
month.loc1<-which(cdata$TIMESTAMP$mday==1&
cdata$TIMESTAMP$hour==0&
cdata$TIMESTAMP$min==0)
month.ticks1 <- seq(cdata$TIMESTAMP[month.loc1[1]],
cdata$TIMESTAMP[month.loc1[length(month.loc1)]],by="months")
png(paste0(plot.path,
"Concord_H_LE_co2_2_",
cdata$TIMESTAMP$year[1]+1900,"_",
cdata$TIMESTAMP$yday[1]+1,"_",
cdata$TIMESTAMP$year[nrow(cdata)]+1900,"_",
cdata$TIMESTAMP$yday[nrow(cdata)]+1,"_",
target.plot.var_1[k1],"_color_",
Sys.Date(),".png"),
width=8,
height=4,
units="in",
res=300,
pointsize = 11,
bg = "white")
par(oma=c(4,4.5,0.5,0.5),mar=c(0,0.5,0,0.5),fig=c(0,0.7,0,1))
plot(cdata$Correct_NR,
cdata[,target.plot.var_1[k1]],
xlab="",
ylab="",
cex=0.5,col=Wd.col_1[WD.SE_LE_co2_by_NR],
las=1,pch=16,
xaxs="i",yaxs="i"
)
mtext(side=2,target.plot.var_1.title[k1],line=3)
mtext(side=1,"Net Radiation",line=2.8)
abline(h=0,col="darkgrey")
axis(1, at = month.ticks1, labels = FALSE, tcl = -0.3)
par(fig=c(0.7,1,0,1),new=T)
hist0A<-hist(cdata[,target.plot.var_1[k1]],
plot=F,nclass=50)
hist1A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==1],
plot=F,breaks=hist0A$breaks)
hist2A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==2],
plot=F,breaks=hist0A$breaks)
barplot(hist1A$counts,
axes=F,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[1],border=NA) # barplot
barplot(hist2A$counts,
axes=F,
add=T,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[2],border=NA) # barplot
legend(0,
length(hist1A$breaks)+1,
fill=Wd.col_1,border=NA,
legend=Wd.legend_1,bty="n",
cex=0.9)
dev.off()
}
#Treatment 120-300, Control 0-120 and 300-360
```r
WD.SE_LE_co2_by_NR<-rep(2,nrow(cdata))
WD.SE_LE_co2_by_NR[which(!is.na(cdata$wind_dir)&
(cdata$wind_dir>120&cdata$wind_dir<=300))]<-1
Wd.legend_1<-c(\120-300 wind \,\0-120 and 300-360 wind\)
Wd.col_1<-c(rgb(0,0,1,0.5,maxColorValue=1),rgb(1,0,0,0.5,maxColorValue=1))
target.plot.var_1<-c(\co2_flux\,\LE\,\H\,\u.\)
target.plot.var_1.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(LE~'('~W~m^{-2}~')'),
expression(H~'('~W~m^{-2}~')'),
expression(u['*']~'('~m~s^{-1}~')'))
for(k1 in 1:length(target.plot.var_1)){
## locate the start of each month
month.loc1<-which(cdata$TIMESTAMP$mday==1&
cdata$TIMESTAMP$hour==0&
cdata$TIMESTAMP$min==0)
month.ticks1 <- seq(cdata$TIMESTAMP[month.loc1[1]],
cdata$TIMESTAMP[month.loc1[length(month.loc1)]],by=\months\)
png(paste0(plot.path,
\Concord_H_LE_co2_3_\,
cdata$TIMESTAMP$year[1]+1900,\_\,
cdata$TIMESTAMP$yday[1]+1,\_\,
cdata$TIMESTAMP$year[nrow(cdata)]+1900,\_\,
cdata$TIMESTAMP$yday[nrow(cdata)]+1,\_\,
target.plot.var_1[k1],\_color_\,
Sys.Date(),\.png\),
width=8,
height=4,
units=\in\,
res=300,
pointsize = 11,
bg = \white\)
par(oma=c(4,4.5,0.5,0.5),mar=c(0,0.5,0,0.5),fig=c(0,0.7,0,1))
plot(cdata$Correct_NR,
cdata[,target.plot.var_1[k1]],
xlab=\\,
ylab=\\,
cex=0.5,col=Wd.col_1[WD.SE_LE_co2_by_NR],
las=1,pch=16,
xaxs=\i\,yaxs=\i\
)
mtext(side=2,target.plot.var_1.title[k1],line=3)
mtext(side=1,\Net Radiation\,line=2.8)
abline(h=0,col=\darkgrey\)
axis(1, at = month.ticks1, labels = FALSE, tcl = -0.3)
par(fig=c(0.7,1,0,1),new=T)
hist0A<-hist(cdata[,target.plot.var_1[k1]],
plot=F,nclass=50)
hist1A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==1],
plot=F,breaks=hist0A$breaks)
hist2A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==2],
plot=F,breaks=hist0A$breaks)
barplot(hist1A$counts,
axes=F,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[1],border=NA) # barplot
barplot(hist2A$counts,
axes=F,
add=T,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[2],border=NA) # barplot
legend(0,
length(hist1A$breaks)+1,
fill=Wd.col_1,border=NA,
legend=Wd.legend_1,bty=\n\,
cex=0.9)
dev.off()
}
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
# Sensible Heat, Latent Heat U* and Co2 fluxes by Net radiation, seperated by wind directions
#Treatment 130-230, Control 230-330
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuV0QuU0VfTEVfY28yX2J5X05SPC1yZXAoMixucm93KGNkYXRhKSlcbldELlNFX0xFX2NvMl9ieV9OUlt3aGljaCghaXMubmEoY2RhdGEkd2luZF9kaXIpJlxuICAgICAgICAgICAgICAgKGNkYXRhJHdpbmRfZGlyPj0xMzAgJiBjZGF0YSR3aW5kX2Rpcjw9MjMwKSYgI3RyZWF0bWVudFxuICAgICAgICAgICAgICAgICFpcy5uYShjZGF0YSR3aW5kX2RpciksXG4gICAgICAgICAgICAgICAgICAgY2RhdGEkd2luZF9kaXI+MjMwICYgY2RhdGEkd2luZF9kaXI8PTMzMCldPC0xICNjb250cm9sXG5XZC5sZWdlbmRfMTwtYyhcXDEzMC0yMzAgd2luZCBcXCxcXDIzMC0zMzAgd2luZFxcKVxuV2QuY29sXzE8LWMocmdiKDAsMCwxLDAuNSxtYXhDb2xvclZhbHVlPTEpLHJnYigxLDAsMCwwLjUsbWF4Q29sb3JWYWx1ZT0xKSlcblxudGFyZ2V0LnBsb3QudmFyXzE8LWMoXFxjbzJfZmx1eFxcLFxcTEVcXCxcXEhcXCxcXHUuXFwpXG50YXJnZXQucGxvdC52YXJfMS50aXRsZTwtYyhleHByZXNzaW9uKEZDficoJ35tdX5tb2x+bV57LTJ9fnNeey0xfX4nKScpLFxuICAgICAgICAgICAgICAgICAgICAgICAgIGV4cHJlc3Npb24oTEV+Jygnfld+bV57LTJ9ficpJyksXG4gICAgICAgICAgICAgICAgICAgICAgICAgZXhwcmVzc2lvbihIficoJ35Xfm1eey0yfX4nKScpLFxuICAgICAgICAgICAgICAgICAgICAgICAgIGV4cHJlc3Npb24odVsnKiddficoJ35tfnNeey0xfX4nKScpKSBcblxuZm9yKGsxIGluIDE6bGVuZ3RoKHRhcmdldC5wbG90LnZhcl8xKSl7XG4gIFxuICAjIyBsb2NhdGUgdGhlIHN0YXJ0IG9mIGVhY2ggbW9udGhcbiAgbW9udGgubG9jMTwtd2hpY2goY2RhdGEkVElNRVNUQU1QJG1kYXk9PTEmXG4gICAgICAgICAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVAkaG91cj09MCZcbiAgICAgICAgICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCRtaW49PTApXG4gIG1vbnRoLnRpY2tzMSA8LSBzZXEoY2RhdGEkVElNRVNUQU1QW21vbnRoLmxvYzFbMV1dLFxuICAgICAgICAgICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QW21vbnRoLmxvYzFbbGVuZ3RoKG1vbnRoLmxvYzEpXV0sYnk9XFxtb250aHNcXClcbiAgXG4gIHBuZyhwYXN0ZTAocGxvdC5wYXRoLFxuICAgICAgICAgICAgIFxcQ29uY29yZF9IX0xFX2NvMl80X1xcLFxuICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCR5ZWFyWzFdKzE5MDAsXFxfXFwsXG4gICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJHlkYXlbMV0rMSxcXF9cXCxcbiAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVAkeWVhcltucm93KGNkYXRhKV0rMTkwMCxcXF9cXCxcbiAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVAkeWRheVtucm93KGNkYXRhKV0rMSxcXF9cXCxcbiAgICAgICAgICAgICB0YXJnZXQucGxvdC52YXJfMVtrMV0sXFxfY29sb3JfXFwsXG4gICAgICAgICAgICAgU3lzLkRhdGUoKSxcXC5wbmdcXCksXG4gICAgICB3aWR0aD04LFxuICAgICAgaGVpZ2h0PTQsXG4gICAgICB1bml0cz1cXGluXFwsXG4gICAgICByZXM9MzAwLFxuICAgICAgcG9pbnRzaXplID0gMTEsXG4gICAgICBiZyA9IFxcd2hpdGVcXClcbiAgXG4gIHBhcihvbWE9Yyg0LDQuNSwwLjUsMC41KSxtYXI9YygwLDAuNSwwLDAuNSksZmlnPWMoMCwwLjcsMCwxKSlcbiAgcGxvdChjZGF0YSRDb3JyZWN0X05SLFxuICAgICAgIGNkYXRhWyx0YXJnZXQucGxvdC52YXJfMVtrMV1dLFxuICAgICAgIHhsYWI9XFxcXCxcbiAgICAgICB5bGFiPVxcXFwsXG4gICAgICAgY2V4PTAuNSxjb2w9V2QuY29sXzFbV0QuU0VfTEVfY28yX2J5X05SXSxcbiAgICAgICBsYXM9MSxwY2g9MTYsXG4gICAgICAgeGF4cz1cXGlcXCx5YXhzPVxcaVxcXG4gIClcbiAgbXRleHQoc2lkZT0yLHRhcmdldC5wbG90LnZhcl8xLnRpdGxlW2sxXSxsaW5lPTMpXG4gIG10ZXh0KHNpZGU9MSxcXE5ldCBSYWRpYXRpb25cXCxsaW5lPTIuOClcbiAgYWJsaW5lKGg9MCxjb2w9XFxkYXJrZ3JleVxcKVxuICBheGlzKDEsIGF0ID0gbW9udGgudGlja3MxLCBsYWJlbHMgPSBGQUxTRSwgdGNsID0gLTAuMylcbiAgXG4gIHBhcihmaWc9YygwLjcsMSwwLDEpLG5ldz1UKVxuICBoaXN0MEE8LWhpc3QoY2RhdGFbLHRhcmdldC5wbG90LnZhcl8xW2sxXV0sXG4gICAgICAgICAgICAgIHBsb3Q9RixuY2xhc3M9NTApXG4gIGhpc3QxQTwtaGlzdChjZGF0YVssdGFyZ2V0LnBsb3QudmFyXzFbazFdXVtXRC5TRV9MRV9jbzJfYnlfTlI9PTFdLFxuICAgICAgICAgICAgICBwbG90PUYsYnJlYWtzPWhpc3QwQSRicmVha3MpXG4gIGhpc3QyQTwtaGlzdChjZGF0YVssdGFyZ2V0LnBsb3QudmFyXzFbazFdXVtXRC5TRV9MRV9jbzJfYnlfTlI9PTJdLFxuICAgICAgICAgICAgICBwbG90PUYsYnJlYWtzPWhpc3QwQSRicmVha3MpIFxuICBcbiAgYmFycGxvdChoaXN0MUEkY291bnRzLFxuICAgICAgICAgIGF4ZXM9RixcbiAgICAgICAgICBob3Jpej1ULFxuICAgICAgICAgIHlsaW09YygwLGxlbmd0aChoaXN0MUEkYnJlYWtzKSsxKSxcbiAgICAgICAgICB4bGltPWMoLTUsbWF4KGMoaGlzdDFBJGNvdW50cyxoaXN0MkEkY291bnRzKSkpLFxuICAgICAgICAgIHNwYWNlPTAsY29sPVdkLmNvbF8xWzFdLGJvcmRlcj1OQSkgIyBiYXJwbG90XG4gIGJhcnBsb3QoaGlzdDJBJGNvdW50cyxcbiAgICAgICAgICBheGVzPUYsXG4gICAgICAgICAgYWRkPVQsXG4gICAgICAgICAgaG9yaXo9VCxcbiAgICAgICAgICB5bGltPWMoMCxsZW5ndGgoaGlzdDFBJGJyZWFrcykrMSksXG4gICAgICAgICAgeGxpbT1jKC01LG1heChjKGhpc3QxQSRjb3VudHMsaGlzdDJBJGNvdW50cykpKSxcbiAgICAgICAgICBzcGFjZT0wLGNvbD1XZC5jb2xfMVsyXSxib3JkZXI9TkEpICMgYmFycGxvdFxuICBsZWdlbmQoMCxcbiAgICAgICAgIGxlbmd0aChoaXN0MUEkYnJlYWtzKSsxLFxuICAgICAgICAgZmlsbD1XZC5jb2xfMSxib3JkZXI9TkEsXG4gICAgICAgICBsZWdlbmQ9V2QubGVnZW5kXzEsYnR5PVxcblxcLFxuICAgICAgICAgY2V4PTAuOSlcbiAgZGV2Lm9mZigpXG59XG4gIFxuYGBgXG5gYGAifQ== -->
```r
```r
WD.SE_LE_co2_by_NR<-rep(2,nrow(cdata))
WD.SE_LE_co2_by_NR[which(!is.na(cdata$wind_dir)&
(cdata$wind_dir>=130 & cdata$wind_dir<=230)& #treatment
!is.na(cdata$wind_dir),
cdata$wind_dir>230 & cdata$wind_dir<=330)]<-1 #control
Wd.legend_1<-c(\130-230 wind \,\230-330 wind\)
Wd.col_1<-c(rgb(0,0,1,0.5,maxColorValue=1),rgb(1,0,0,0.5,maxColorValue=1))
target.plot.var_1<-c(\co2_flux\,\LE\,\H\,\u.\)
target.plot.var_1.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(LE~'('~W~m^{-2}~')'),
expression(H~'('~W~m^{-2}~')'),
expression(u['*']~'('~m~s^{-1}~')'))
for(k1 in 1:length(target.plot.var_1)){
## locate the start of each month
month.loc1<-which(cdata$TIMESTAMP$mday==1&
cdata$TIMESTAMP$hour==0&
cdata$TIMESTAMP$min==0)
month.ticks1 <- seq(cdata$TIMESTAMP[month.loc1[1]],
cdata$TIMESTAMP[month.loc1[length(month.loc1)]],by=\months\)
png(paste0(plot.path,
\Concord_H_LE_co2_4_\,
cdata$TIMESTAMP$year[1]+1900,\_\,
cdata$TIMESTAMP$yday[1]+1,\_\,
cdata$TIMESTAMP$year[nrow(cdata)]+1900,\_\,
cdata$TIMESTAMP$yday[nrow(cdata)]+1,\_\,
target.plot.var_1[k1],\_color_\,
Sys.Date(),\.png\),
width=8,
height=4,
units=\in\,
res=300,
pointsize = 11,
bg = \white\)
par(oma=c(4,4.5,0.5,0.5),mar=c(0,0.5,0,0.5),fig=c(0,0.7,0,1))
plot(cdata$Correct_NR,
cdata[,target.plot.var_1[k1]],
xlab=\\,
ylab=\\,
cex=0.5,col=Wd.col_1[WD.SE_LE_co2_by_NR],
las=1,pch=16,
xaxs=\i\,yaxs=\i\
)
mtext(side=2,target.plot.var_1.title[k1],line=3)
mtext(side=1,\Net Radiation\,line=2.8)
abline(h=0,col=\darkgrey\)
axis(1, at = month.ticks1, labels = FALSE, tcl = -0.3)
par(fig=c(0.7,1,0,1),new=T)
hist0A<-hist(cdata[,target.plot.var_1[k1]],
plot=F,nclass=50)
hist1A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==1],
plot=F,breaks=hist0A$breaks)
hist2A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==2],
plot=F,breaks=hist0A$breaks)
barplot(hist1A$counts,
axes=F,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[1],border=NA) # barplot
barplot(hist2A$counts,
axes=F,
add=T,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[2],border=NA) # barplot
legend(0,
length(hist1A$breaks)+1,
fill=Wd.col_1,border=NA,
legend=Wd.legend_1,bty=\n\,
cex=0.9)
dev.off()
}
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
# Sensible Heat, Latent Heat U* and Co2 fluxes by Net radiation, seperated by wind directions
#Treatment 150-225, Control 230-300
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuV0QuU0VfTEVfY28yX2J5X05SPC1yZXAoMixucm93KGNkYXRhKSlcbldELlNFX0xFX2NvMl9ieV9OUlt3aGljaCghaXMubmEoY2RhdGEkd2luZF9kaXIpJlxuICAgICAgICAgICAgICAgKGNkYXRhJHdpbmRfZGlyPj0xNTAgJiBjZGF0YSR3aW5kX2Rpcjw9MjI1KSYgI3RyZWF0bWVudFxuICAgICAgICAgICAgICAgICFpcy5uYShjZGF0YSR3aW5kX2RpciksXG4gICAgICAgICAgICAgICAgICAgY2RhdGEkd2luZF9kaXI+PTIzMCAmIGNkYXRhJHdpbmRfZGlyPD0zMDApXTwtMSAjY29udHJvbFxuV2QubGVnZW5kXzE8LWMoXCIxNTAtMjI1IHdpbmQgXCIsXCIyMzAtMzAwIHdpbmRcIilcbldkLmNvbF8xPC1jKHJnYigwLDAsMSwwLjUsbWF4Q29sb3JWYWx1ZT0xKSxyZ2IoMSwwLDAsMC41LG1heENvbG9yVmFsdWU9MSkpXG5cbnRhcmdldC5wbG90LnZhcl8xPC1jKFwiY28yX2ZsdXhcIixcIkxFXCIsXCJIXCIsXCJ1LlwiKVxudGFyZ2V0LnBsb3QudmFyXzEudGl0bGU8LWMoZXhwcmVzc2lvbihGQ34nKCd+bXV+bW9sfm1eey0yfX5zXnstMX1+JyknKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICBleHByZXNzaW9uKExFficoJ35Xfm1eey0yfX4nKScpLFxuICAgICAgICAgICAgICAgICAgICAgICAgIGV4cHJlc3Npb24oSH4nKCd+V35tXnstMn1+JyknKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICBleHByZXNzaW9uKHVbJyonXX4nKCd+bX5zXnstMX1+JyknKSkgXG5cbmZvcihrMSBpbiAxOmxlbmd0aCh0YXJnZXQucGxvdC52YXJfMSkpe1xuICBcbiAgIyMgbG9jYXRlIHRoZSBzdGFydCBvZiBlYWNoIG1vbnRoXG4gIG1vbnRoLmxvYzE8LXdoaWNoKGNkYXRhJFRJTUVTVEFNUCRtZGF5PT0xJlxuICAgICAgICAgICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJGhvdXI9PTAmXG4gICAgICAgICAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVAkbWluPT0wKVxuICBtb250aC50aWNrczEgPC0gc2VxKGNkYXRhJFRJTUVTVEFNUFttb250aC5sb2MxWzFdXSxcbiAgICAgICAgICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUFttb250aC5sb2MxW2xlbmd0aChtb250aC5sb2MxKV1dLGJ5PVwibW9udGhzXCIpXG4gIFxuICBwbmcocGFzdGUwKHBsb3QucGF0aCxcbiAgICAgICAgICAgICBcIkNvbmNvcmRfSF9MRV9jbzJfNV9cIixcbiAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVAkeWVhclsxXSsxOTAwLFwiX1wiLFxuICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCR5ZGF5WzFdKzEsXCJfXCIsXG4gICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJHllYXJbbnJvdyhjZGF0YSldKzE5MDAsXCJfXCIsXG4gICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJHlkYXlbbnJvdyhjZGF0YSldKzEsXCJfXCIsXG4gICAgICAgICAgICAgdGFyZ2V0LnBsb3QudmFyXzFbazFdLFwiX2NvbG9yX1wiLFxuICAgICAgICAgICAgIFN5cy5EYXRlKCksXCIucG5nXCIpLFxuICAgICAgd2lkdGg9OCxcbiAgICAgIGhlaWdodD00LFxuICAgICAgdW5pdHM9XCJpblwiLFxuICAgICAgcmVzPTMwMCxcbiAgICAgIHBvaW50c2l6ZSA9IDExLFxuICAgICAgYmcgPSBcIndoaXRlXCIpXG4gIFxuICBwYXIob21hPWMoNCw0LjUsMC41LDAuNSksbWFyPWMoMCwwLjUsMCwwLjUpLGZpZz1jKDAsMC43LDAsMSkpXG4gIHBsb3QoY2RhdGEkQ29ycmVjdF9OUixcbiAgICAgICBjZGF0YVssdGFyZ2V0LnBsb3QudmFyXzFbazFdXSxcbiAgICAgICB4bGFiPVwiXCIsXG4gICAgICAgeWxhYj1cIlwiLFxuICAgICAgIGNleD0wLjUsY29sPVdkLmNvbF8xW1dELlNFX0xFX2NvMl9ieV9OUl0sXG4gICAgICAgbGFzPTEscGNoPTE2LFxuICAgICAgIHhheHM9XCJpXCIseWF4cz1cImlcIlxuICApXG4gIG10ZXh0KHNpZGU9Mix0YXJnZXQucGxvdC52YXJfMS50aXRsZVtrMV0sbGluZT0zKVxuICBtdGV4dChzaWRlPTEsXCJOZXQgUmFkaWF0aW9uXCIsbGluZT0yLjgpXG4gIGFibGluZShoPTAsY29sPVwiZGFya2dyZXlcIilcbiAgYXhpcygxLCBhdCA9IG1vbnRoLnRpY2tzMSwgbGFiZWxzID0gRkFMU0UsIHRjbCA9IC0wLjMpXG4gIFxuICBwYXIoZmlnPWMoMC43LDEsMCwxKSxuZXc9VClcbiAgaGlzdDBBPC1oaXN0KGNkYXRhWyx0YXJnZXQucGxvdC52YXJfMVtrMV1dLFxuICAgICAgICAgICAgICBwbG90PUYsbmNsYXNzPTUwKVxuICBoaXN0MUE8LWhpc3QoY2RhdGFbLHRhcmdldC5wbG90LnZhcl8xW2sxXV1bV0QuU0VfTEVfY28yX2J5X05SPT0xXSxcbiAgICAgICAgICAgICAgcGxvdD1GLGJyZWFrcz1oaXN0MEEkYnJlYWtzKVxuICBoaXN0MkE8LWhpc3QoY2RhdGFbLHRhcmdldC5wbG90LnZhcl8xW2sxXV1bV0QuU0VfTEVfY28yX2J5X05SPT0yXSxcbiAgICAgICAgICAgICAgcGxvdD1GLGJyZWFrcz1oaXN0MEEkYnJlYWtzKSBcbiAgXG4gIGJhcnBsb3QoaGlzdDFBJGNvdW50cyxcbiAgICAgICAgICBheGVzPUYsXG4gICAgICAgICAgaG9yaXo9VCxcbiAgICAgICAgICB5bGltPWMoMCxsZW5ndGgoaGlzdDFBJGJyZWFrcykrMSksXG4gICAgICAgICAgeGxpbT1jKC01LG1heChjKGhpc3QxQSRjb3VudHMsaGlzdDJBJGNvdW50cykpKSxcbiAgICAgICAgICBzcGFjZT0wLGNvbD1XZC5jb2xfMVsxXSxib3JkZXI9TkEpICMgYmFycGxvdFxuICBiYXJwbG90KGhpc3QyQSRjb3VudHMsXG4gICAgICAgICAgYXhlcz1GLFxuICAgICAgICAgIGFkZD1ULFxuICAgICAgICAgIGhvcml6PVQsXG4gICAgICAgICAgeWxpbT1jKDAsbGVuZ3RoKGhpc3QxQSRicmVha3MpKzEpLFxuICAgICAgICAgIHhsaW09YygtNSxtYXgoYyhoaXN0MUEkY291bnRzLGhpc3QyQSRjb3VudHMpKSksXG4gICAgICAgICAgc3BhY2U9MCxjb2w9V2QuY29sXzFbMl0sYm9yZGVyPU5BKSAjIGJhcnBsb3RcbiAgbGVnZW5kKDAsXG4gICAgICAgICBsZW5ndGgoaGlzdDFBJGJyZWFrcykrMSxcbiAgICAgICAgIGZpbGw9V2QuY29sXzEsYm9yZGVyPU5BLFxuICAgICAgICAgbGVnZW5kPVdkLmxlZ2VuZF8xLGJ0eT1cIm5cIixcbiAgICAgICAgIGNleD0wLjkpXG4gIGRldi5vZmYoKVxufVxuICBcbmBgYCJ9 -->
```r
WD.SE_LE_co2_by_NR<-rep(2,nrow(cdata))
WD.SE_LE_co2_by_NR[which(!is.na(cdata$wind_dir)&
(cdata$wind_dir>=150 & cdata$wind_dir<=225)& #treatment
!is.na(cdata$wind_dir),
cdata$wind_dir>=230 & cdata$wind_dir<=300)]<-1 #control
Wd.legend_1<-c("150-225 wind ","230-300 wind")
Wd.col_1<-c(rgb(0,0,1,0.5,maxColorValue=1),rgb(1,0,0,0.5,maxColorValue=1))
target.plot.var_1<-c("co2_flux","LE","H","u.")
target.plot.var_1.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'),
expression(LE~'('~W~m^{-2}~')'),
expression(H~'('~W~m^{-2}~')'),
expression(u['*']~'('~m~s^{-1}~')'))
for(k1 in 1:length(target.plot.var_1)){
## locate the start of each month
month.loc1<-which(cdata$TIMESTAMP$mday==1&
cdata$TIMESTAMP$hour==0&
cdata$TIMESTAMP$min==0)
month.ticks1 <- seq(cdata$TIMESTAMP[month.loc1[1]],
cdata$TIMESTAMP[month.loc1[length(month.loc1)]],by="months")
png(paste0(plot.path,
"Concord_H_LE_co2_5_",
cdata$TIMESTAMP$year[1]+1900,"_",
cdata$TIMESTAMP$yday[1]+1,"_",
cdata$TIMESTAMP$year[nrow(cdata)]+1900,"_",
cdata$TIMESTAMP$yday[nrow(cdata)]+1,"_",
target.plot.var_1[k1],"_color_",
Sys.Date(),".png"),
width=8,
height=4,
units="in",
res=300,
pointsize = 11,
bg = "white")
par(oma=c(4,4.5,0.5,0.5),mar=c(0,0.5,0,0.5),fig=c(0,0.7,0,1))
plot(cdata$Correct_NR,
cdata[,target.plot.var_1[k1]],
xlab="",
ylab="",
cex=0.5,col=Wd.col_1[WD.SE_LE_co2_by_NR],
las=1,pch=16,
xaxs="i",yaxs="i"
)
mtext(side=2,target.plot.var_1.title[k1],line=3)
mtext(side=1,"Net Radiation",line=2.8)
abline(h=0,col="darkgrey")
axis(1, at = month.ticks1, labels = FALSE, tcl = -0.3)
par(fig=c(0.7,1,0,1),new=T)
hist0A<-hist(cdata[,target.plot.var_1[k1]],
plot=F,nclass=50)
hist1A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==1],
plot=F,breaks=hist0A$breaks)
hist2A<-hist(cdata[,target.plot.var_1[k1]][WD.SE_LE_co2_by_NR==2],
plot=F,breaks=hist0A$breaks)
barplot(hist1A$counts,
axes=F,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[1],border=NA) # barplot
barplot(hist2A$counts,
axes=F,
add=T,
horiz=T,
ylim=c(0,length(hist1A$breaks)+1),
xlim=c(-5,max(c(hist1A$counts,hist2A$counts))),
space=0,col=Wd.col_1[2],border=NA) # barplot
legend(0,
length(hist1A$breaks)+1,
fill=Wd.col_1,border=NA,
legend=Wd.legend_1,bty="n",
cex=0.9)
dev.off()
}
#treatment 150-225, Control 230-300. My way of overlaying just to make sure above code is legit #does Net radiation (categorical independent variable) have an impact on co2 fluxes #covariate is treatment v.s control wind directions
```r
plot(
cdata$Correct_NR[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$daytime == 0 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$daytime == 0 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-20, 20),
xlim = c(0, 600),
xlab = 'Net Radiation',
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
main = 'Night time fluxes by Net Radiation',
cex = 0.6,
col = \blue\
)
par(new = TRUE)
plot(
cdata$Correct_NR[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$daytime == 0 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$daytime == 0 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-20, 20),
xlim = c(0, 600),
xlab = 'Net Radiation',
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
#main='Night time fluxes by Net Radiation',
cex = 0.6,
col = \red\
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#all time fluxes by Net radiation filtering based on u*> 0.1
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuXG5wbG90KFxuICBjZGF0YSRDb3JyZWN0X05SW2NkYXRhJHdpbmRfZGlyID49IDE1MCAmXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAyMjUgJlxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICBcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGNkYXRhJHUuID4gMC4xXSxcbiAgY2RhdGEkY28yX2ZsdXhbY2RhdGEkd2luZF9kaXIgPj0gMTUwICZcbiAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAyMjUgICYgY2RhdGEkdS4gPiAwLjFdLFxuICAjIHhheHQ9ICduJyxcbiAgIyB5YXh0PSAnbicsXG4gIGFibGluZShoID0gMCwgY29sID0gXCJkYXJrZ3JleVwiKSxcbiAgeWxpbSA9IGMoLTIwLCAyMCksXG4gIHhsaW0gPSBjKC0xMDAsIDcwMCksXG4gIHhsYWIgPSAnTmV0IFJhZGlhdGlvbicsXG4gIHlsYWIgPSBleHByZXNzaW9uKEZDIH4gJygnIH4gbXUgfiBtb2wgfiBtIF4ge1xuICAgIC0yXG4gIH0gfiBzIF4ge1xuICAgIC0xXG4gIH0gfiAnKScpLFxuICBjZXgubGFiPTAuOCxcbiAgXG4gIG1haW4gPSAnQWxsIGZsdXhlcyBieSBOZXQgUmFkaWF0aW9uJyxcbiAgXG4gIGNleCA9IDAuNixcbiAgY29sID0gXCJibHVlXCJcbilcblxucGFyKG5ldyA9IFRSVUUpXG5gYGAifQ== -->
```r
plot(
cdata$Correct_NR[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(-100, 700),
xlab = 'Net Radiation',
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
cex.lab=0.8,
main = 'All fluxes by Net Radiation',
cex = 0.6,
col = "blue"
)
par(new = TRUE)
plot(
cdata$Correct_NR[cdata$wind_dir >= 230 &
cdata$wind_dir <= 330 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 230 &
cdata$wind_dir <= 330 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(-100, 700),
cex.lab= 0.8,
xlab = 'Net Radiation',
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
#main='Night time fluxes by Net Radiation',
cex = 0.6,
col = "red"
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
#ANCOVA comparison of co2 fluxes by NR for the treatment and control areas
#With ustar filter of greater than 0.1
#ANCOVA of different wind directions
#Day and Night sorting for ANCOVA analysis
cdata$day_and_night_wind_filter<-rep("exclude",nrow(cdata))
cdata$day_and_night_wind_filter[which(!is.na(cdata$wind_dir)&
cdata$wind_dir>=150&
cdata$wind_dir<=225&
cdata$u. > 0.1)]<-"treatment"
cdata$day_and_night_wind_filter[which(!is.na(cdata$wind_dir)&
cdata$wind_dir>=230&
cdata$wind_dir<=300&
cdata$u. > 0.1)]<-"control"
table(cdata$day_and_night_wind_filter)
control exclude treatment
7796 21657 15785
##############################################################################
#nightime sorting of data for ANCOVA analysis
cdata$night_wind_filter<-rep("exclude",nrow(cdata))
cdata$night_wind_filter[which(!is.na(cdata$wind_dir)&
cdata$wind_dir>=150&
cdata$wind_dir<=225 &
cdata$u. > 0.1 &
cdata$daytime== 0)]<-"treatment"
cdata$night_wind_filter[which(!is.na(cdata$wind_dir)&
cdata$wind_dir>=230&
cdata$wind_dir<=300 &
cdata$u. > 0.1&
cdata$daytime== 0)]<-"control"
table(cdata$night_wind_filter)
control exclude treatment
2690 31835 10713
###########################################################################################
#daytime sorting of data for ANCOVA analsi
cdata$day_wind_filter<-rep("exclude",nrow(cdata))
cdata$day_wind_filter[which(!is.na(cdata$wind_dir)&
cdata$wind_dir>=150&
cdata$wind_dir<=225 &
cdata$u. > 0.1 &
cdata$daytime== 1)]<-"treatment"
cdata$day_wind_filter[which(!is.na(cdata$wind_dir)&
cdata$wind_dir>=230&
cdata$wind_dir<=300 &
cdata$u. > 0.1&
cdata$daytime== 1)]<-"control"
table(cdata$day_wind_filter)
control exclude treatment
5106 35060 5072
####################################################################################
#ANCOVA of LE by NR>0
#No significant difference
LE_by_wind_NR <-lmer(LE ~ Correct_NR + day_and_night_wind_filter +(1|day_and_night_wind_filter),
data = cdata[cdata$Correct_NR>0&
cdata$day_and_night_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),])
summary(LE_by_wind_NR)
Linear mixed model fit by REML ['lmerMod']
Formula:
LE ~ Correct_NR + day_and_night_wind_filter + (1 | day_and_night_wind_filter)
Data:
cdata[cdata$Correct_NR > 0 & cdata$day_and_night_wind_filter !=
"exclude" & (cdata$TIMESTAMP$mon <= 3 | cdata$TIMESTAMP$mon >=
11), ]
REML criterion at convergence: 20153.3
Scaled residuals:
Min 1Q Median 3Q Max
-2.2811 -0.7413 -0.0890 0.6243 4.6670
Random effects:
Groups Name Variance Std.Dev.
day_and_night_wind_filter (Intercept) 25.21 5.021
Residual 1350.14 36.744
Number of obs: 2006, groups: day_and_night_wind_filter, 2
Fixed effects:
Estimate Std. Error
(Intercept) 43.566759 5.305525
Correct_NR 0.109566 0.005267
day_and_night_wind_filtertreatment -2.840260 7.300572
t value
(Intercept) 8.212
Correct_NR 20.802
day_and_night_wind_filtertreatment -0.389
Correlation of Fixed Effects:
(Intr) Crr_NR
Correct_NR -0.248
dy_nd_ngh__ -0.695 0.055
Anova(LE_by_wind_NR)
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: LE
Chisq Df Pr(>Chisq)
Correct_NR 432.7089 1 <2e-16 ***
day_and_night_wind_filter 0.1514 1 0.6972
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#ANCOVA of LE by NR neg and pos (day and night time)
#Don't do ANCOVA for night time radiation. just noise
# summary(lm(LE ~ Correct_NR + day_and_night_wind_filter,
# data = cdata[cdata$Correct_NR&
# cdata$day_and_night_wind_filter&cdata$DOY!="exclude"&
# (cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),]))
#ANCOVA of co2 flux by NR with NR>0 (day time)
#No significant difference
co2_by_NR_wind_dir <- lmer(co2_flux ~ Correct_NR + day_and_night_wind_filter +(1|day_and_night_wind_filter),
data = cdata[cdata$Correct_NR>0&
cdata$day_and_night_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),])
Warning: unable to evaluate scaled gradientWarning: Hessian is numerically singular: parameters are not uniquely determined
summary(co2_by_NR_wind_dir)
Linear mixed model fit by REML ['lmerMod']
Formula:
co2_flux ~ Correct_NR + day_and_night_wind_filter + (1 | day_and_night_wind_filter)
Data:
cdata[cdata$Correct_NR > 0 & cdata$day_and_night_wind_filter !=
"exclude" & (cdata$TIMESTAMP$mon <= 3 | cdata$TIMESTAMP$mon >=
11), ]
REML criterion at convergence: 11928.7
Scaled residuals:
Min 1Q Median 3Q Max
-6.9959 -0.4417 0.1891 0.5588 5.5212
Random effects:
Groups Name Variance Std.Dev.
day_and_night_wind_filter (Intercept) 0.03625 0.1904
Residual 25.41332 5.0412
Number of obs: 1962, groups: day_and_night_wind_filter, 2
Fixed effects:
Estimate Std. Error
(Intercept) -2.545271 0.305241
Correct_NR -0.006992 0.000728
day_and_night_wind_filtertreatment -0.909966 0.357843
t value
(Intercept) -8.339
Correct_NR -9.604
day_and_night_wind_filtertreatment -2.543
Correlation of Fixed Effects:
(Intr) Crr_NR
Correct_NR -0.603
dy_nd_ngh__ -0.633 0.151
optimizer (nloptwrap) convergence code: 0 (OK)
unable to evaluate scaled gradient
Hessian is numerically singular: parameters are not uniquely determined
Anova(co2_by_NR_wind_dir)
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: co2_flux
Chisq Df Pr(>Chisq)
Correct_NR 92.2436 1 < 2e-16 ***
day_and_night_wind_filter 6.4664 1 0.01099 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
shapiro.test(resid(co2_by_NR_wind_dir))
Shapiro-Wilk normality test
data: resid(co2_by_NR_wind_dir)
W = 0.92168, p-value < 2.2e-16
#ANCOVA of co2 flux by NR with NR neg and pos (day and night time)
#No significant difference
# summary(lm(co2_flux ~ Correct_NR + day_and_night_wind_filter,
# data = cdata[cdata$Correct_NR&
# cdata$day_and_night_wind_filter!="exclude"&
# (cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),]))
#ANCOVA of H by NR
# No Significant differance accounting for dependence between treatment and control sites
H_by_wind_dir_NR <-lmer(H ~ Correct_NR + day_and_night_wind_filter +(1|day_and_night_wind_filter),
data = cdata[cdata$Correct_NR>0&
cdata$day_and_night_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),])
summary(H_by_wind_dir_NR)
Linear mixed model fit by REML ['lmerMod']
Formula:
H ~ Correct_NR + day_and_night_wind_filter + (1 | day_and_night_wind_filter)
Data:
cdata[cdata$Correct_NR > 0 & cdata$day_and_night_wind_filter !=
"exclude" & (cdata$TIMESTAMP$mon <= 3 | cdata$TIMESTAMP$mon >=
11), ]
REML criterion at convergence: 24492
Scaled residuals:
Min 1Q Median 3Q Max
-3.4524 -0.6980 0.0011 0.6417 3.8754
Random effects:
Groups Name Variance Std.Dev.
day_and_night_wind_filter (Intercept) 177.4 13.32
Residual 2199.8 46.90
Number of obs: 2325, groups: day_and_night_wind_filter, 2
Fixed effects:
Estimate Std. Error
(Intercept) -24.63321 13.47537
Correct_NR 0.50566 0.00617
day_and_night_wind_filtertreatment 11.02477 18.94309
t value
(Intercept) -1.828
Correct_NR 81.961
day_and_night_wind_filtertreatment 0.582
Correlation of Fixed Effects:
(Intr) Crr_NR
Correct_NR -0.117
dy_nd_ngh__ -0.705 0.027
Anova(H_by_wind_dir_NR)
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: H
Chisq Df Pr(>Chisq)
Correct_NR 6717.5696 1 <2e-16 ***
day_and_night_wind_filter 0.3387 1 0.5606
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#ANCOVA of u* by NR
#Noi significant difference
ustar_by_NR_wind_dir <- lmer(u. ~ Correct_NR + day_and_night_wind_filter + (1|day_and_night_wind_filter),
data = cdata[cdata$Correct_NR>0&
cdata$day_and_night_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),])
summary(ustar_by_NR_wind_dir)
Linear mixed model fit by REML ['lmerMod']
Formula:
u. ~ Correct_NR + day_and_night_wind_filter + (1 | day_and_night_wind_filter)
Data:
cdata[cdata$Correct_NR > 0 & cdata$day_and_night_wind_filter !=
"exclude" & (cdata$TIMESTAMP$mon <= 3 | cdata$TIMESTAMP$mon >=
11), ]
REML criterion at convergence: -4114.1
Scaled residuals:
Min 1Q Median 3Q Max
-2.3513 -0.6913 -0.1251 0.5533 4.1694
Random effects:
Groups Name Variance Std.Dev.
day_and_night_wind_filter (Intercept) 0.001704 0.04128
Residual 0.010443 0.10219
Number of obs: 2408, groups: day_and_night_wind_filter, 2
Fixed effects:
Estimate Std. Error
(Intercept) 1.950e-01 4.151e-02
Correct_NR 2.222e-04 1.321e-05
day_and_night_wind_filtertreatment 1.057e-01 5.853e-02
t value
(Intercept) 4.699
Correct_NR 16.817
day_and_night_wind_filtertreatment 1.805
Correlation of Fixed Effects:
(Intr) Crr_NR
Correct_NR -0.080
dy_nd_ngh__ -0.706 0.018
Anova(ustar_by_NR_wind_dir)
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: u.
Chisq Df Pr(>Chisq)
Correct_NR 282.8054 1 < 2e-16 ***
day_and_night_wind_filter 3.2581 1 0.07107 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
##########################################################################
#ANCOVA Nightime fluxes by air temperature
nighttime_co2_by_airT <- lmer(co2_flux ~ air_temperature_adj + night_wind_filter +(1|night_wind_filter),
data = cdata[cdata$air_temperature_adj&
cdata$night_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),])
summary(nighttime_co2_by_airT)
Linear mixed model fit by REML ['lmerMod']
Formula:
co2_flux ~ air_temperature_adj + night_wind_filter + (1 | night_wind_filter)
Data:
cdata[cdata$air_temperature_adj & cdata$night_wind_filter !=
"exclude" & (cdata$TIMESTAMP$mon <= 3 | cdata$TIMESTAMP$mon >=
11), ]
REML criterion at convergence: 20043.1
Scaled residuals:
Min 1Q Median 3Q Max
-17.5512 -0.4429 -0.1114 0.3362 5.0929
Random effects:
Groups Name Variance Std.Dev.
night_wind_filter (Intercept) 0.01563 0.125
Residual 10.13489 3.184
Number of obs: 3887, groups: night_wind_filter, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 2.89774 0.22227 13.037
air_temperature_adj -0.03080 0.01366 -2.255
night_wind_filtertreatment 0.29814 0.22389 1.332
Correlation of Fixed Effects:
(Intr) ar_tm_
ar_tmprtr_d -0.609
nght_wnd_fl -0.668 0.071
Anova(nighttime_co2_by_airT)
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: co2_flux
Chisq Df Pr(>Chisq)
air_temperature_adj 5.0865 1 0.02411 *
night_wind_filter 1.7733 1 0.18298
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#ANCOVA daytime fluxes by air temperature
day_co2_by_airT <- lmer(co2_flux ~ air_temperature_adj + day_wind_filter +(1|day_wind_filter),
data = cdata[cdata$air_temperature_adj&
cdata$day_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),])
Warning: unable to evaluate scaled gradientWarning: Hessian is numerically singular: parameters are not uniquely determined
summary(day_co2_by_airT)
Linear mixed model fit by REML ['lmerMod']
Formula:
co2_flux ~ air_temperature_adj + day_wind_filter + (1 | day_wind_filter)
Data:
cdata[cdata$air_temperature_adj & cdata$day_wind_filter != "exclude" &
(cdata$TIMESTAMP$mon <= 3 | cdata$TIMESTAMP$mon >= 11), ]
REML criterion at convergence: 14310.3
Scaled residuals:
Min 1Q Median 3Q Max
-6.6430 -0.4716 0.1159 0.6260 5.1184
Random effects:
Groups Name Variance Std.Dev.
day_wind_filter (Intercept) 0.5493 0.7412
Residual 29.8933 5.4675
Number of obs: 2294, groups: day_wind_filter, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) -2.26713 0.85058 -2.665
air_temperature_adj -0.15163 0.02606 -5.819
day_wind_filtertreatment 0.94276 1.07707 0.875
Correlation of Fixed Effects:
(Intr) ar_tm_
ar_tmprtr_d -0.449
dy_wnd_fltr -0.670 0.088
optimizer (nloptwrap) convergence code: 0 (OK)
unable to evaluate scaled gradient
Hessian is numerically singular: parameters are not uniquely determined
Anova(day_co2_by_airT)
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: co2_flux
Chisq Df Pr(>Chisq)
air_temperature_adj 33.8651 1 5.907e-09 ***
day_wind_filter 0.7662 1 0.3814
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#ANCOVA day and night fluxes by air temperature. Stupid to do. because interdependences between temperature and night and just being respiration temp drops at night where we will just see respiration
###################################################################################
#ANCOVA Nightime fluxes by soil temperature
night_co2_by_soil_t1 <-lmer(co2_flux ~ TC_Avg.1. + night_wind_filter + (1|night_wind_filter),
data = cdata[cdata$TC_Avg.1.&
cdata$night_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),])
Warning: unable to evaluate scaled gradientWarning: Hessian is numerically singular: parameters are not uniquely determined
summary(night_co2_by_soil_t1)
Linear mixed model fit by REML ['lmerMod']
Formula:
co2_flux ~ TC_Avg.1. + night_wind_filter + (1 | night_wind_filter)
Data:
cdata[cdata$TC_Avg.1. & cdata$night_wind_filter != "exclude" &
(cdata$TIMESTAMP$mon <= 3 | cdata$TIMESTAMP$mon >= 11), ]
REML criterion at convergence: 19999.6
Scaled residuals:
Min 1Q Median 3Q Max
-17.7396 -0.4341 -0.1116 0.3303 5.0625
Random effects:
Groups Name Variance Std.Dev.
night_wind_filter (Intercept) 0.01546 0.1243
Residual 10.02114 3.1656
Number of obs: 3887, groups: night_wind_filter, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 3.60385 0.22699 15.877
TC_Avg.1. -0.07817 0.01114 -7.017
night_wind_filtertreatment 0.28716 0.22216 1.293
Correlation of Fixed Effects:
(Intr) TC_A.1
TC_Avg.1. -0.635
nght_wnd_fl -0.629 0.030
optimizer (nloptwrap) convergence code: 0 (OK)
unable to evaluate scaled gradient
Hessian is numerically singular: parameters are not uniquely determined
Anova(night_co2_by_soil_t1)
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: co2_flux
Chisq Df Pr(>Chisq)
TC_Avg.1. 49.2341 1 2.272e-12 ***
night_wind_filter 1.6707 1 0.1962
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#ANCOVA daytime fluxes by soil temperature
day_co2_by_soil_t1 <-lmer(co2_flux ~ TC_Avg.1. + day_wind_filter + (1|day_wind_filter),
data = cdata[cdata$TC_Avg.1.&
cdata$day_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),])
summary(day_co2_by_soil_t1)
Linear mixed model fit by REML ['lmerMod']
Formula:
co2_flux ~ TC_Avg.1. + day_wind_filter + (1 | day_wind_filter)
Data:
cdata[cdata$TC_Avg.1. & cdata$day_wind_filter != "exclude" &
(cdata$TIMESTAMP$mon <= 3 | cdata$TIMESTAMP$mon >= 11), ]
REML criterion at convergence: 14332.3
Scaled residuals:
Min 1Q Median 3Q Max
-6.6214 -0.4806 0.0750 0.6228 4.9226
Random effects:
Groups Name Variance Std.Dev.
day_wind_filter (Intercept) 0.3718 0.6097
Residual 30.3396 5.5081
Number of obs: 2292, groups: day_wind_filter, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) -4.581846 0.715875 -6.400
TC_Avg.1. 0.005373 0.020460 0.263
day_wind_filtertreatment 1.521293 0.895116 1.700
Correlation of Fixed Effects:
(Intr) TC_A.1
TC_Avg.1. -0.468
dy_wnd_fltr -0.660 0.074
Anova(day_co2_by_soil_t1)
Analysis of Deviance Table (Type II Wald chisquare tests)
Response: co2_flux
Chisq Df Pr(>Chisq)
TC_Avg.1. 0.0690 1 0.79285
day_wind_filter 2.8885 1 0.08922 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
##############################################################################################
#ANCOVA Nightime fluxes by soil VWC
summary(lm(co2_flux ~ VWC_Avg + night_wind_filter,
data = cdata[cdata$VWC_Avg&
cdata$night_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),]))
Call:
lm(formula = co2_flux ~ VWC_Avg + night_wind_filter, data = cdata[cdata$VWC_Avg &
cdata$night_wind_filter != "exclude" & (cdata$TIMESTAMP$mon <=
3 | cdata$TIMESTAMP$mon >= 11), ])
Residuals:
Min 1Q Median 3Q Max
-55.409 -1.289 -0.319 0.990 15.262
Coefficients:
Estimate Std. Error t value
(Intercept) 1.8119 0.1405 12.89
VWC_Avg 3.4992 0.2793 12.53
night_wind_filtertreatment 0.2902 0.1369 2.12
Pr(>|t|)
(Intercept) <2e-16 ***
VWC_Avg <2e-16 ***
night_wind_filtertreatment 0.0341 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 3.116 on 3699 degrees of freedom
(1472 observations deleted due to missingness)
Multiple R-squared: 0.04244, Adjusted R-squared: 0.04192
F-statistic: 81.96 on 2 and 3699 DF, p-value: < 2.2e-16
#ANCOVA daytime fluxes by soil temperature
summary(lm(co2_flux ~ VWC_Avg + day_wind_filter,
data = cdata[cdata$VWC_Avg&
cdata$day_wind_filter!="exclude"&
(cdata$TIMESTAMP$mon<=3|cdata$TIMESTAMP$mon>=11),]))
Call:
lm(formula = co2_flux ~ VWC_Avg + day_wind_filter, data = cdata[cdata$VWC_Avg &
cdata$day_wind_filter != "exclude" & (cdata$TIMESTAMP$mon <=
3 | cdata$TIMESTAMP$mon >= 11), ])
Residuals:
Min 1Q Median 3Q Max
-35.772 -2.781 0.129 3.254 26.347
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -3.2509 0.2181 -14.902 < 2e-16
VWC_Avg -6.2116 0.6353 -9.777 < 2e-16
day_wind_filtertreatment 1.6152 0.2274 7.103 1.64e-12
(Intercept) ***
VWC_Avg ***
day_wind_filtertreatment ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 5.374 on 2245 degrees of freedom
(720 observations deleted due to missingness)
Multiple R-squared: 0.06014, Adjusted R-squared: 0.0593
F-statistic: 71.83 on 2 and 2245 DF, p-value: < 2.2e-16
#LE by Net Radiation
#expression(LE~'('~W~m^{-2}~')')
```r
plot(
cdata$Correct_NR[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$u. > 0.1],
cdata$LE[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-50, 300),
xlim = c(-100, 700),
xlab = 'Net Radiation',
ylab = expression(LE~'('~W~m^{-2}~')'),
main = 'Latent Heat by Net Radiation',
cex = 0.6,
col = \blue\
)
par(new = TRUE)
plot(
cdata$Correct_NR[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$u. > 0.1],
cdata$LE[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-50, 300),
xlim = c(-100, 700),
xlab = 'Net Radiation',
ylab =expression(LE~'('~W~m^{-2}~')'),
#main='Night time fluxes by Net Radiation',
cex = 0.6,
col = \red\
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#sensible heat by net radiation
#expression(H~'('~W~m^{-2}~')')
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuXG5wbG90KFxuICBjZGF0YSRDb3JyZWN0X05SW2NkYXRhJHdpbmRfZGlyID49IDE1MCAmXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAyMjUgJlxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICBcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGNkYXRhJHUuID4gMC4xXSxcbiAgY2RhdGEkSFtjZGF0YSR3aW5kX2RpciA+PSAxNTAgJlxuICAgICAgICAgICAgICAgICAgIGNkYXRhJHdpbmRfZGlyIDw9IDIyNSAgJiBjZGF0YSR1LiA+IDAuMV0sXG4gICMgeGF4dD0gJ24nLFxuICAjIHlheHQ9ICduJyxcbiAgYWJsaW5lKGggPSAwLCBjb2wgPSBcXGRhcmtncmV5XFwpLFxuICB5bGltID0gYygtMTAwLCA0MDApLFxuICB4bGltID0gYygtMTAwLCA3MDApLFxuICB4bGFiID0gJ05ldCBSYWRpYXRpb24nLFxuICB5bGFiID0gZXhwcmVzc2lvbihIficoJ35Xfm1eey0yfX4nKScpLFxuICBcbiAgbWFpbiA9ICdTZW5zaWJsZSBIZWF0IGJ5IE5ldCBSYWRpYXRpb24nLFxuICBcbiAgY2V4ID0gMC42LFxuICBjb2wgPSBcXGJsdWVcXFxuKVxuXG5wYXIobmV3ID0gVFJVRSlcblxuXG5wbG90KFxuICBjZGF0YSRDb3JyZWN0X05SW2NkYXRhJHdpbmRfZGlyID49IDIzMCAmXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAzMDAgJlxuICAgICAgICAgICAgICAgICAgICAgICAgICAgIFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY2RhdGEkdS4gPiAwLjFdLFxuICBjZGF0YSRIW2NkYXRhJHdpbmRfZGlyID49IDIzMCAmXG4gICAgICAgICAgICAgICAgICAgY2RhdGEkd2luZF9kaXIgPD0gMzAwICYgY2RhdGEkdS4gPiAwLjFdLFxuICAjIHhheHQ9ICduJyxcbiAgIyB5YXh0PSAnbicsXG4gIGFibGluZShoID0gMCwgY29sID0gXFxkYXJrZ3JleVxcKSxcbiAgeWxpbSA9IGMoLTEwMCwgNDAwKSxcbiAgeGxpbSA9IGMoLTEwMCwgNzAwKSxcbiAgeGxhYiA9ICdOZXQgUmFkaWF0aW9uJyxcbiAgeWxhYiA9ZXhwcmVzc2lvbihIficoJ35Xfm1eey0yfX4nKScpLFxuICBcbiAgXG4gIFxuICBjZXggPSAwLjYsXG4gIGNvbCA9IFxccmVkXFxcbilcblxuI2NyZWF0aW5nIHRoZSBsZWdlbmRcbmxlZ2VuZChcbiAgeCA9ICdib3R0b21yaWdodCcsXG4gIGxlZ2VuZCA9IGMoJ1RyZWF0bWVudCAoMTUwLTIyNSknLCAnQ29udHJvbCgyMzAtMzAwKScpLFxuICBwY2ggPSAxLFxuICBjb2wgPSBjKCdibHVlJywgJ3JlZCcpXG4pXG5cbmBgYFxuYGBgIn0= -->
```r
```r
plot(
cdata$Correct_NR[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$u. > 0.1],
cdata$H[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-100, 400),
xlim = c(-100, 700),
xlab = 'Net Radiation',
ylab = expression(H~'('~W~m^{-2}~')'),
main = 'Sensible Heat by Net Radiation',
cex = 0.6,
col = \blue\
)
par(new = TRUE)
plot(
cdata$Correct_NR[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$u. > 0.1],
cdata$H[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-100, 400),
xlim = c(-100, 700),
xlab = 'Net Radiation',
ylab =expression(H~'('~W~m^{-2}~')'),
cex = 0.6,
col = \red\
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#U* by net radiation
#expression(u['*']~'('~m~s^{-1}~')')
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxucGxvdChcbiAgY2RhdGEkQ29ycmVjdF9OUltjZGF0YSR3aW5kX2RpciA+PSAxNTAgJlxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY2RhdGEkd2luZF9kaXIgPD0gMjI1ICZcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjZGF0YSR1LiA+IDAuMV0sXG4gIGNkYXRhJHUuW2NkYXRhJHdpbmRfZGlyID49IDE1MCAmXG4gICAgICAgICAgICAgICAgICAgY2RhdGEkd2luZF9kaXIgPD0gMjI1ICAmIGNkYXRhJHUuID4gMC4xXSxcbiAgIyB4YXh0PSAnbicsXG4gICMgeWF4dD0gJ24nLFxuICBhYmxpbmUoaCA9IDAsIGNvbCA9IFxcZGFya2dyZXlcXCksXG4gIHlsaW0gPSBjKDAuMSwgMSksXG4gIHhsaW0gPSBjKC0xMDAsIDcwMCksXG4gIHhsYWIgPSAnTmV0IFJhZGlhdGlvbicsXG4gIHlsYWIgPSBleHByZXNzaW9uKHVbJyonXX4nKCd+bX5zXnstMX1+JyknKSxcbiAgXG4gIG1haW4gPSAnVSogYnkgTmV0IFJhZGlhdGlvbicsXG4gIFxuICBjZXggPSAwLjYsXG4gIGNvbCA9IFxcYmx1ZVxcXG4pXG5cbnBhcihuZXcgPSBUUlVFKVxuXG5cbnBsb3QoXG4gIGNkYXRhJENvcnJlY3RfTlJbY2RhdGEkd2luZF9kaXIgPj0gMjMwICZcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGNkYXRhJHdpbmRfZGlyIDw9IDMwMCAmXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjZGF0YSR1LiA+IDAuMV0sXG4gIGNkYXRhJHUuW2NkYXRhJHdpbmRfZGlyID49IDIzMCAmXG4gICAgICAgICAgICAgICAgICAgY2RhdGEkd2luZF9kaXIgPD0gMzAwICYgY2RhdGEkdS4gPiAwLjFdLFxuICAjIHhheHQ9ICduJyxcbiAgIyB5YXh0PSAnbicsXG4gIGFibGluZShoID0gMCwgY29sID0gXFxkYXJrZ3JleVxcKSxcbiAgeWxpbSA9IGMoMC4xLCAxKSxcbiAgeGxpbSA9IGMoLTEwMCwgNzAwKSxcbiAgeGxhYiA9ICdOZXQgUmFkaWF0aW9uJyxcbiAgeWxhYiA9ZXhwcmVzc2lvbih1WycqJ11+Jygnfm1+c157LTF9ficpJyksXG4gIFxuICBcbiAgXG4gIGNleCA9IDAuNixcbiAgY29sID0gXFxyZWRcXFxuKVxuXG4jY3JlYXRpbmcgdGhlIGxlZ2VuZFxubGVnZW5kKFxuICB4ID0gJ2JvdHRvbXJpZ2h0JyxcbiAgbGVnZW5kID0gYygnVHJlYXRtZW50ICgxNTAtMjI1KScsICdDb250cm9sKDIzMC0zMDApJyksXG4gIHBjaCA9IDEsXG4gIGNvbCA9IGMoJ2JsdWUnLCAncmVkJylcbilcbmBgYFxuYGBgIn0= -->
```r
```r
plot(
cdata$Correct_NR[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$u. > 0.1],
cdata$u.[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(0.1, 1),
xlim = c(-100, 700),
xlab = 'Net Radiation',
ylab = expression(u['*']~'('~m~s^{-1}~')'),
main = 'U* by Net Radiation',
cex = 0.6,
col = \blue\
)
par(new = TRUE)
plot(
cdata$Correct_NR[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$u. > 0.1],
cdata$u.[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(0.1, 1),
xlim = c(-100, 700),
xlab = 'Net Radiation',
ylab =expression(u['*']~'('~m~s^{-1}~')'),
cex = 0.6,
col = \red\
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#looking at break down of fluxes by wind direction
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuXG5cbndpbmRfMTwtaGlzdChjZGF0YSR3aW5kX2RpcltjZGF0YSR3aW5kX2Rpcj49IDE1MCAmXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAyMjUgJiBjZGF0YSR1LiA+IDAuMV0sIFxuICAgICBicmVha3MgPSA3LCBcbiAgICAgeGxpbSA9IGMoMCwzNjApLCBcbiAgICAgeWxpbSA9IGMoMCwyNjAwKSwgXG4gICAgIGNvbCA9ICdibHVlJyxcbiAgICAgeGxhYiA9ICAnV2luZCBEaXJlY3Rpb24nLFxuICAgICAjbGFiZWxzID0gVFJVRSxcbiAgICAgbWFpbiA9IFwiV2luZCBEaXJlY3Rpb24gRnJlcXVlbmN5XCIpXG5cblxucGFyKG5ldz1UUlVFKVxuYGBgIn0= -->
```r
wind_1<-hist(cdata$wind_dir[cdata$wind_dir>= 150 &
cdata$wind_dir <= 225 & cdata$u. > 0.1],
breaks = 7,
xlim = c(0,360),
ylim = c(0,2600),
col = 'blue',
xlab = 'Wind Direction',
#labels = TRUE,
main = "Wind Direction Frequency")
par(new=TRUE)
hist(cdata$wind_dir[cdata$wind_dir>= 230 &
cdata$wind_dir <= 300 & cdata$u. > 0.1],
breaks = 7,
xlim = c(0,360),
ylim = c(0,2600),
col = 'red',
xlab = "",
#labels = TRUE,
main = "")
par(new=TRUE)
hist(cdata$wind_dir[cdata$wind_dir>= 0 &
cdata$wind_dir <150 & cdata$u. > 0.1],
breaks = 15,
xlim = c(0,360),
ylim = c(0,2600),
xlab = "",
#labels = TRUE,
main = "")
par(new=TRUE)
hist(cdata$wind_dir[cdata$wind_dir> 300 & cdata$u. > 0.1 ],
breaks = 6,
xlim = c(0,360),
ylim = c(0,2600),
xlab = "",
#labels = TRUE,
main="")
par(new=TRUE)
hist(cdata$wind_dir[cdata$wind_dir> 225 &
cdata$wind_dir<230 & cdata$u. > 0.1],
breaks = 1,
xlim = c(0,360),
ylim = c(0,2600),
col = "green",
xlab = "",
#labels = TRUE,
main="")
legend(
x = 'topright',
legend = c('Treatment (150-225)', 'Control(230-300)', 'gap(225-230)', 'excluded'),
lty = 1,
col = c('blue', 'red', 'green', 'grey')
)
#treatment data points
freq_table(cdata$wind_dir[cdata$wind_dir>= 150 &
cdata$wind_dir <= 225 & cdata$u. > 0.1])
freq_table(cdata$wind_dir[cdata$wind_dir>= 230 &
cdata$wind_dir <= 300 & cdata$u. > 0.1])
#excluded data (includes wetlands, etc.)
freq_table(cdata$wind_dir[cdata$wind_dir<150 |
cdata$wind_dir > 300 & cdata$u. > 0.1])
#excluded dat to create a buffer between treatment and contrl
freq_table(cdata$wind_dir[cdata$wind_dir> 225 &
cdata$wind_dir < 230 & cdata$u. > 0.1])
# wind_table<-freq_table(cdata$wind_dir[cdata$wind_dir>= 150 &
# cdata$wind_dir <= 225]&
# cdata$wind_dir[cdata$wind_dir>= 230 &
# cdata$wind_dir <= 300])
#
# wind_table
#night time co2 fluxes by air temperature # daytime==’0 #using wind directions for treatment 150-225, Control 230-300 #im struggling with this.
```r
# WD.co2_night_by_air_temp<-rep(NA,nrow(cdata)) #keep just in case you have missing values
# WD.co2_night_by_air_temp[which(!is.na(cdata$wind_dir)&
# (cdata$wind_dir>=150 & cdata$wind_dir<=225)&
# !is.na(cdata$daytime)&
# cdata$daytime == '0')]<-1 #night time treatment
#
# WD.co2_night_by_air_temp[which(!is.na(cdata$wind_dir)&
# (cdata$wind_dir>=230 | cdata$wind_dir<=300)&
# !is.na(cdata$daytime)&
# cdata$daytime == '0')]<-2 #night time control
#
# WD.co2_night_by_air_temp[which(!is.na(cdata$wind_dir)&
# (cdata$wind_dir>270&cdata$wind_dir<=360)&
# !is.na(cdata$daytime)&
# cdata$daytime == '1')]<-3 #daytime time treatment
#
# WD.co2_night_by_air_temp[which(!is.na(cdata$wind_dir)&
# (cdata$wind_dir<270 | cdata$wind_dir>=360)&
# !is.na(cdata$daytime)&
# cdata$daytime == '1')]<-4 #day time control
#
# table(WD.co2_night_by_air_temp)
#
#
# #could create a night time cdata
# #WD.co2_night_by_air_temp[which(!is.na(cdata$daytime)&
# # (cdata$daytime == '0' & cdata$daytime=='1'))]<-2 #changed to 2 just now
#
# WD.legend<-c(\Treatment wind\,\Control wind\)
# WD.col<-c(rgb(0,0,1,0.5,maxColorValue=1),rgb(1,0,0,0.5,maxColorValue=1)) #four groups have to assign four colors. can assign whiteto hide what we don't see.
#
# target.plot.var<-c(\co2_flux\)
# target.plot.var.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'))
#
# for(k1 in 1:length(target.plot.var)){
#
# ## locate the start of each month
# month.loc<-which(cdata$TIMESTAMP$mday==1&
# cdata$TIMESTAMP$hour==0&
# cdata$TIMESTAMP$min==0)
# month.ticks <- seq(cdata$TIMESTAMP[month.loc[1]],
# cdata$TIMESTAMP[month.loc[length(month.loc)]],by=\months\)
#
# png(paste0(plot.path,
# \Concord_8_\,
# cdata$TIMESTAMP$year[1]+1900,\_\,
# cdata$TIMESTAMP$yday[1]+1,\_\,
# cdata$TIMESTAMP$year[nrow(cdata)]+1900,\_\,
# cdata$TIMESTAMP$yday[nrow(cdata)]+1,\_\,
# target.plot.var[k1],\_color_\,
# Sys.Date(),\.png\),
# width=8,
# height=4,
# units=\in\,
# res=300,
# pointsize = 11,
# bg = \white\)
#
# par(oma=c(4,4.5,0.5,0.5),mar=c(0,0.5,0,0.5),fig=c(0,0.7,0,1))
# plot(cdata$air_temperature_adj[], #plotting everything, just making different color. Can hide group 4& 3
# cdata[,target.plot.var[k1]],#add square bracket here.
# xlab=\\,
# ylab=\\,
# cex=0.5,col=WD.col[1],
# las=1,pch=16,
# xaxs=\i\,yaxs=\i\
#
# )
#
# #use a new function called point Cdata:air temperature and target group, pch. plotting and adding anoter level. have two groups so just need 1
# mtext(side=2,target.plot.var.title[k1],line=3)
# mtext(side=1,\Air Temperature Night\,line=2.8)
# abline(h=0,col=\darkgrey\)
# #axis(1, at = month.ticks, labels = FALSE, tcl = -0.3)
#
# par(fig=c(0.7,1,0,1),new=T)
# hist0<-hist(cdata[,target.plot.var[k1]],
# plot=F,nclass=50)
# hist1<-hist(cdata[,target.plot.var[k1]][WD.grp==1], #subsetting and only using half of the data
# plot=F,breaks=hist0$breaks)
# hist2<-hist(cdata[,target.plot.var[k1]][WD.grp==2],
# plot=F,breaks=hist0$breaks)
#
# barplot(hist1$counts,
# axes=F,
# horiz=T,
# ylim=c(0,length(hist1$breaks)+1),
# xlim=c(-5,max(c(hist1$counts,hist2$counts))),
# space=0,col=WD.col[1],border=NA) # barplot
# barplot(hist2$counts,
# axes=F,
# add=T,
# horiz=T,
# ylim=c(0,length(hist1$breaks)+1),
# xlim=c(-5,max(c(hist1$counts,hist2$counts))),
# space=0,col=WD.col[2],border=NA) # barplot
# legend(0,
# length(hist1$breaks)+1,
# fill=WD.col,border=NA,
# legend=WD.legend,bty=\n\,
# cex=0.9)
# dev.off()
# }
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#co2 fluxes by night time. using my plotting over method, PAR new=TRUE
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuXG5wbG90KFxuICBjZGF0YSRhaXJfdGVtcGVyYXR1cmVfYWRqW2NkYXRhJHdpbmRfZGlyID49IDE1MCAmXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAyMjUgJlxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY2RhdGEkZGF5dGltZSA9PSAwICZcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGNkYXRhJHUuID4gMC4xXSxcbiAgY2RhdGEkY28yX2ZsdXhbY2RhdGEkd2luZF9kaXIgPj0gMTUwICZcbiAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAyMjUgJiBjZGF0YSRkYXl0aW1lID09IDAgJiBjZGF0YSR1LiA+IDAuMV0sXG4gICMgeGF4dD0gJ24nLFxuICAjIHlheHQ9ICduJyxcbiAgYWJsaW5lKGggPSAwLCBjb2wgPSBcImRhcmtncmV5XCIpLFxuICB5bGltID0gYygtMjAsIDIwKSxcbiAgeGxpbSA9IGMoMCwgMzUpLFxuICB4bGFiID0gZXhwcmVzc2lvbihBaXIgfiB0ZW1wZXJhdHVyZSB+ICcoJyB+IGRlZ3JlZSB+IEMgfiAnKScpLFxuICB5bGFiID0gZXhwcmVzc2lvbihGQyB+ICcoJyB+IG11IH4gbW9sIH4gbSBeIHtcbiAgICAtMlxuICB9IH4gcyBeIHtcbiAgICAtMVxuICB9IH4gJyknKSxcbiAgXG4gIG1haW4gPSAnTmlnaHQgdGltZSBmbHV4ZXMgYnkgYWlyIHRlbXBlcmF0dXJlJyxcbiAgXG4gIGNleCA9IDAuNixcbiAgY29sID0gXCJibHVlXCJcbilcblxucGFyKG5ldyA9IFRSVUUpXG5gYGAifQ== -->
```r
plot(
cdata$air_temperature_adj[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$daytime == 0 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$daytime == 0 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(0, 35),
xlab = expression(Air ~ temperature ~ '(' ~ degree ~ C ~ ')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
main = 'Night time fluxes by air temperature',
cex = 0.6,
col = "blue"
)
par(new = TRUE)
plot(
cdata$air_temperature_adj[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$daytime == 0 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$daytime == 0 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(0, 35),
xlab = expression(Air ~ temperature ~ '(' ~ degree ~ C ~ ')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
#main='Night time fluxes by air temperature',
cex = 0.6,
col = "red"
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
NA
NA
NA
#daytime co2 fluxes
plot(
cdata$air_temperature_adj[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$daytime == 1 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$daytime == 1 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(0, 35),
xlab = expression(Air ~ temperature ~ '(' ~ degree ~ C ~ ')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
main = 'Day time fluxes by air temperature',
cex = 0.6,
col = "blue"
)
par(new = TRUE)
plot(
cdata$air_temperature_adj[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$daytime == 1 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$daytime == 1 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(0, 35),
xlab = expression(Air ~ temperature ~ '(' ~ degree ~ C ~ ')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
#main='Day time fluxes by air temperature',
cex = 0.6,
col = "red"
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
#nighttime fluxes by soil temperature
plot(
cdata$TC_Avg.1.[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$daytime == 0 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$daytime == 0 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(0, 35),
xlab = expression(Soil ~ Surface ~ temperature ~ '(' ~ degree ~ C ~
')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
main = 'Night time fluxes by surface soil temperature (TC1)',
cex = 0.6,
col = "blue"
)
par(new = TRUE)
plot(
cdata$TC_Avg.1.[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$daytime == 0 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$daytime == 0 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(0, 35),
xlab = expression(Soil ~ Surface ~ temperature ~ '(' ~ degree ~ C ~
')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
#main='Night time fluxes by Surface Soil Temperature',
cex = 0.6,
col = "red"
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
#Daytime fluxes by soil temperature TC1
plot(
cdata$TC_Avg.1.[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$daytime == 1 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$daytime == 1 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(0, 35),
xlab = expression(Soil ~ Surface ~ temperature ~ '(' ~ degree ~ C ~
')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
main = 'Day time fluxes by surface soil temperature (TC1)',
cex = 0.6,
col = "blue"
)
par(new = TRUE)
plot(
cdata$TC_Avg.1.[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$daytime == 1 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$daytime == 1 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = "darkgrey"),
ylim = c(-20, 20),
xlim = c(0, 35),
xlab = expression(Soil ~ Surface ~ temperature ~ '(' ~ degree ~ C ~
')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
#main='Night time fluxes by Surface Soil Temperature',
cex = 0.6,
col = "red"
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
#night time fluxes by VWc (soil moisture)
```r
plot(
cdata$VWC_Avg[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$daytime == 0 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$daytime == 0 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-20, 20),
xlim = c(0, 1),
xlab = expression(VWC ~ '(' ~ Percentage ~ ')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
main = 'Night time fluxes by VWC',
cex = 0.6,
col = \blue\
)
par(new = TRUE)
plot(
cdata$VWC_Avg[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$daytime == 0 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$daytime == 0 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-20, 20),
xlim = c(0, 1),
xlab = expression(VWC ~ '(' ~ Percentage ~ ')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
#main='Night time fluxes by Surface Soil Temperature',
cex = 0.6,
col = \red\
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#daytime fluxes by VWC
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuXG5wbG90KFxuICBjZGF0YSRWV0NfQXZnW2NkYXRhJHdpbmRfZGlyID49IDE1MCAmXG4gICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAyMjUgJlxuICAgICAgICAgICAgICAgICAgY2RhdGEkZGF5dGltZSA9PSAxICZcbiAgICAgICAgICAgICAgICAgIGNkYXRhJHUuID4gMC4xXSxcbiAgY2RhdGEkY28yX2ZsdXhbY2RhdGEkd2luZF9kaXIgPj0gMTUwICZcbiAgICAgICAgICAgICAgICAgICBjZGF0YSR3aW5kX2RpciA8PSAyMjUgJiBjZGF0YSRkYXl0aW1lID09IDEgJiBjZGF0YSR1LiA+IDAuMV0sXG4gICMgeGF4dD0gJ24nLFxuICAjIHlheHQ9ICduJyxcbiAgYWJsaW5lKGggPSAwLCBjb2wgPSBcXGRhcmtncmV5XFwpLFxuICB5bGltID0gYygtMjAsIDIwKSxcbiAgeGxpbSA9IGMoMCwgMSksXG4gIHhsYWIgPSBleHByZXNzaW9uKFZXQyB+ICcoJyB+IFBlcmNlbnRhZ2UgfiAnKScpLFxuICB5bGFiID0gZXhwcmVzc2lvbihGQyB+ICcoJyB+IG11IH4gbW9sIH4gbSBeIHtcbiAgICAtMlxuICB9IH4gcyBeIHtcbiAgICAtMVxuICB9IH4gJyknKSxcbiAgXG4gIG1haW4gPSAnRGF5IHRpbWUgZmx1eGVzIGJ5IFZXQycsXG4gIFxuICBjZXggPSAwLjYsXG4gIGNvbCA9IFxcYmx1ZVxcXG4pXG5cbnBhcihuZXcgPSBUUlVFKVxuXG5cbnBsb3QoXG4gIGNkYXRhJFZXQ19BdmdbY2RhdGEkd2luZF9kaXIgPj0gMjMwICZcbiAgICAgICAgICAgICAgICAgIGNkYXRhJHdpbmRfZGlyIDw9IDMwMCAmXG4gICAgICAgICAgICAgICAgICBjZGF0YSRkYXl0aW1lID09IDEgJlxuICAgICAgICAgICAgICAgICAgY2RhdGEkdS4gPiAwLjFdLFxuICBjZGF0YSRjbzJfZmx1eFtjZGF0YSR3aW5kX2RpciA+PSAyMzAgJlxuICAgICAgICAgICAgICAgICAgIGNkYXRhJHdpbmRfZGlyIDw9IDMwMCAmIGNkYXRhJGRheXRpbWUgPT0gMSAmIGNkYXRhJHUuID4gMC4xXSxcbiAgIyB4YXh0PSAnbicsXG4gICMgeWF4dD0gJ24nLFxuICBhYmxpbmUoaCA9IDAsIGNvbCA9IFxcZGFya2dyZXlcXCksXG4gIHlsaW0gPSBjKC0yMCwgMjApLFxuICB4bGltID0gYygwLCAxKSxcbiAgeGxhYiA9IGV4cHJlc3Npb24oVldDIH4gJygnIH4gUGVyY2VudGFnZSB+ICcpJyksXG4gIHlsYWIgPSBleHByZXNzaW9uKEZDIH4gJygnIH4gbXUgfiBtb2wgfiBtIF4ge1xuICAgIC0yXG4gIH0gfiBzIF4ge1xuICAgIC0xXG4gIH0gfiAnKScpLFxuICBcbiAgI21haW49J0RheSB0aW1lIGZsdXhlcyBieSBTdXJmYWNlIFNvaWwgVGVtcGVyYXR1cmUnLFxuICBcbiAgY2V4ID0gMC42LFxuICBjb2wgPSBcXHJlZFxcXG4pXG5cbiNjcmVhdGluZyB0aGUgbGVnZW5kXG5sZWdlbmQoXG4gIHggPSAnYm90dG9tcmlnaHQnLFxuICBsZWdlbmQgPSBjKCdUcmVhdG1lbnQgKDE1MC0yMjUpJywgJ0NvbnRyb2woMjMwLTMwMCknKSxcbiAgcGNoID0gMSxcbiAgY29sID0gYygnYmx1ZScsICdyZWQnKVxuKVxuXG5gYGBcbmBgYCJ9 -->
```r
```r
plot(
cdata$VWC_Avg[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 &
cdata$daytime == 1 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 150 &
cdata$wind_dir <= 225 & cdata$daytime == 1 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-20, 20),
xlim = c(0, 1),
xlab = expression(VWC ~ '(' ~ Percentage ~ ')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
main = 'Day time fluxes by VWC',
cex = 0.6,
col = \blue\
)
par(new = TRUE)
plot(
cdata$VWC_Avg[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 &
cdata$daytime == 1 &
cdata$u. > 0.1],
cdata$co2_flux[cdata$wind_dir >= 230 &
cdata$wind_dir <= 300 & cdata$daytime == 1 & cdata$u. > 0.1],
# xaxt= 'n',
# yaxt= 'n',
abline(h = 0, col = \darkgrey\),
ylim = c(-20, 20),
xlim = c(0, 1),
xlab = expression(VWC ~ '(' ~ Percentage ~ ')'),
ylab = expression(FC ~ '(' ~ mu ~ mol ~ m ^ {
-2
} ~ s ^ {
-1
} ~ ')'),
#main='Day time fluxes by Surface Soil Temperature',
cex = 0.6,
col = \red\
)
#creating the legend
legend(
x = 'bottomright',
legend = c('Treatment (150-225)', 'Control(230-300)'),
pch = 1,
col = c('blue', 'red')
)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#co2 fluxes by air temp, night and day
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuV0QuU0VfTEVfY28yX2J5X05SPC1yZXAoMixucm93KGNkYXRhKSlcbldELlNFX0xFX2NvMl9ieV9OUlt3aGljaCghaXMubmEoY2RhdGEkd2luZF9kaXIpJlxuICAgICAgICAgICAgICAgKGNkYXRhJHdpbmRfZGlyPj0xNTAgJiBjZGF0YSR3aW5kX2Rpcjw9MjI1KSYgI3RyZWF0bWVudFxuICAgICAgICAgICAgICAgICFpcy5uYShjZGF0YSR3aW5kX2RpciksXG4gICAgICAgICAgICAgICAgICAgY2RhdGEkd2luZF9kaXI+PTIzMCAmIGNkYXRhJHdpbmRfZGlyPD0zMDApXTwtMSAjY29udHJvbFxuV2QubGVnZW5kXzE8LWMoXFwxNTAtMjI1IHdpbmQgXFwsXFwyMzAtMzAwIHdpbmRcXClcbldELmNvbDwtYyhyZ2IoMCwwLDEsMC41LG1heENvbG9yVmFsdWU9MSkscmdiKDEsMCwwLDAuNSxtYXhDb2xvclZhbHVlPTEpKVxuXG50YXJnZXQucGxvdC52YXI8LWMoXFxjbzJfZmx1eFxcKVxudGFyZ2V0LnBsb3QudmFyLnRpdGxlPC1jKGV4cHJlc3Npb24oRkN+Jygnfm11fm1vbH5tXnstMn1+c157LTF9ficpJykpIFxuXG5mb3IoazEgaW4gMTpsZW5ndGgodGFyZ2V0LnBsb3QudmFyKSl7XG4gIFxuICAjIyBsb2NhdGUgdGhlIHN0YXJ0IG9mIGVhY2ggbW9udGhcbiAgbW9udGgubG9jPC13aGljaChjZGF0YSRUSU1FU1RBTVAkbWRheT09MSZcbiAgICAgICAgICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCRob3VyPT0wJlxuICAgICAgICAgICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJG1pbj09MClcbiAgbW9udGgudGlja3MgPC0gc2VxKGNkYXRhJFRJTUVTVEFNUFttb250aC5sb2NbMV1dLFxuICAgICAgICAgICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QW21vbnRoLmxvY1tsZW5ndGgobW9udGgubG9jKV1dLGJ5PVxcbW9udGhzXFwpXG4gIFxuICBwbmcocGFzdGUwKHBsb3QucGF0aCxcbiAgICAgICAgICAgICBcXENvbmNvcmRfOV9cXCxcbiAgICAgICAgICAgICBjZGF0YSRUSU1FU1RBTVAkeWVhclsxXSsxOTAwLFxcX1xcLFxuICAgICAgICAgICAgIGNkYXRhJFRJTUVTVEFNUCR5ZGF5WzFdKzEsXFxfXFwsXG4gICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJHllYXJbbnJvdyhjZGF0YSldKzE5MDAsXFxfXFwsXG4gICAgICAgICAgICAgY2RhdGEkVElNRVNUQU1QJHlkYXlbbnJvdyhjZGF0YSldKzEsXFxfXFwsXG4gICAgICAgICAgICAgdGFyZ2V0LnBsb3QudmFyW2sxXSxcXF9jb2xvcl9cXCxcbiAgICAgICAgICAgICBTeXMuRGF0ZSgpLFxcLnBuZ1xcKSxcbiAgICAgIHdpZHRoPTgsXG4gICAgICBoZWlnaHQ9NCxcbiAgICAgIHVuaXRzPVxcaW5cXCxcbiAgICAgIHJlcz0zMDAsXG4gICAgICBwb2ludHNpemUgPSAxMSxcbiAgICAgIGJnID0gXFx3aGl0ZVxcKVxuICBcbiAgcGFyKG9tYT1jKDQsNC41LDAuNSwwLjUpLG1hcj1jKDAsMC41LDAsMC41KSxmaWc9YygwLDAuNywwLDEpKVxuICBwbG90KGNkYXRhJGFpcl90ZW1wZXJhdHVyZV9hZGosXG4gICAgICAgY2RhdGFbLHRhcmdldC5wbG90LnZhcltrMV1dLFxuICAgICAgIHhsYWI9XFxcXCxcbiAgICAgICB5bGFiPVxcXFwsXG4gICAgICAgY2V4PTAuNSxjb2w9V0QuY29sW1dELmdycF0sXG4gICAgICAgbGFzPTEscGNoPTE2LFxuICAgICAgIHhheHM9XFxpXFwseWF4cz1cXGlcXFxuICApXG4gIG10ZXh0KHNpZGU9Mix0YXJnZXQucGxvdC52YXIudGl0bGVbazFdLGxpbmU9MylcbiAgbXRleHQoc2lkZT0xLFxcQWlyIFRlbXBlcmF0dXJlIERheSBhbmQgTmlnaHRcXCxsaW5lPTIuOClcbiAgYWJsaW5lKGg9MCxjb2w9XFxkYXJrZ3JleVxcKVxuICAjYXhpcygxLCBhdCA9IG1vbnRoLnRpY2tzLCBsYWJlbHMgPSBGQUxTRSwgdGNsID0gLTAuMylcbiAgXG4gIHBhcihmaWc9YygwLjcsMSwwLDEpLG5ldz1UKVxuICBoaXN0MDwtaGlzdChjZGF0YVssdGFyZ2V0LnBsb3QudmFyW2sxXV0sXG4gICAgICAgICAgICAgIHBsb3Q9RixuY2xhc3M9NTApXG4gIGhpc3QxPC1oaXN0KGNkYXRhWyx0YXJnZXQucGxvdC52YXJbazFdXVtXRC5ncnA9PTFdLFxuICAgICAgICAgICAgICBwbG90PUYsYnJlYWtzPWhpc3QwJGJyZWFrcylcbiAgaGlzdDI8LWhpc3QoY2RhdGFbLHRhcmdldC5wbG90LnZhcltrMV1dW1dELmdycD09Ml0sXG4gICAgICAgICAgICAgIHBsb3Q9RixicmVha3M9aGlzdDAkYnJlYWtzKSBcbiAgXG4gIGJhcnBsb3QoaGlzdDEkY291bnRzLFxuICAgICAgICAgIGF4ZXM9RixcbiAgICAgICAgICBob3Jpej1ULFxuICAgICAgICAgIHlsaW09YygwLGxlbmd0aChoaXN0MSRicmVha3MpKzEpLFxuICAgICAgICAgIHhsaW09YygtNSxtYXgoYyhoaXN0MSRjb3VudHMsaGlzdDIkY291bnRzKSkpLFxuICAgICAgICAgIHNwYWNlPTAsY29sPVdELmNvbFsxXSxib3JkZXI9TkEpICMgYmFycGxvdFxuICBiYXJwbG90KGhpc3QyJGNvdW50cyxcbiAgICAgICAgICBheGVzPUYsXG4gICAgICAgICAgYWRkPVQsXG4gICAgICAgICAgaG9yaXo9VCxcbiAgICAgICAgICB5bGltPWMoMCxsZW5ndGgoaGlzdDEkYnJlYWtzKSsxKSxcbiAgICAgICAgICB4bGltPWMoLTUsbWF4KGMoaGlzdDEkY291bnRzLGhpc3QyJGNvdW50cykpKSxcbiAgICAgICAgICBzcGFjZT0wLGNvbD1XRC5jb2xbMl0sYm9yZGVyPU5BKSAjIGJhcnBsb3RcbiAgbGVnZW5kKDAsXG4gICAgICAgICBsZW5ndGgoaGlzdDEkYnJlYWtzKSsxLFxuICAgICAgICAgZmlsbD1XRC5jb2wsYm9yZGVyPU5BLFxuICAgICAgICAgbGVnZW5kPVdELmxlZ2VuZCxidHk9XFxuXFwsXG4gICAgICAgICBjZXg9MC45KVxuICBkZXYub2ZmKClcbn1cbmBgYFxuYGBgIn0= -->
```r
```r
WD.SE_LE_co2_by_NR<-rep(2,nrow(cdata))
WD.SE_LE_co2_by_NR[which(!is.na(cdata$wind_dir)&
(cdata$wind_dir>=150 & cdata$wind_dir<=225)& #treatment
!is.na(cdata$wind_dir),
cdata$wind_dir>=230 & cdata$wind_dir<=300)]<-1 #control
Wd.legend_1<-c(\150-225 wind \,\230-300 wind\)
WD.col<-c(rgb(0,0,1,0.5,maxColorValue=1),rgb(1,0,0,0.5,maxColorValue=1))
target.plot.var<-c(\co2_flux\)
target.plot.var.title<-c(expression(FC~'('~mu~mol~m^{-2}~s^{-1}~')'))
for(k1 in 1:length(target.plot.var)){
## locate the start of each month
month.loc<-which(cdata$TIMESTAMP$mday==1&
cdata$TIMESTAMP$hour==0&
cdata$TIMESTAMP$min==0)
month.ticks <- seq(cdata$TIMESTAMP[month.loc[1]],
cdata$TIMESTAMP[month.loc[length(month.loc)]],by=\months\)
png(paste0(plot.path,
\Concord_9_\,
cdata$TIMESTAMP$year[1]+1900,\_\,
cdata$TIMESTAMP$yday[1]+1,\_\,
cdata$TIMESTAMP$year[nrow(cdata)]+1900,\_\,
cdata$TIMESTAMP$yday[nrow(cdata)]+1,\_\,
target.plot.var[k1],\_color_\,
Sys.Date(),\.png\),
width=8,
height=4,
units=\in\,
res=300,
pointsize = 11,
bg = \white\)
par(oma=c(4,4.5,0.5,0.5),mar=c(0,0.5,0,0.5),fig=c(0,0.7,0,1))
plot(cdata$air_temperature_adj,
cdata[,target.plot.var[k1]],
xlab=\\,
ylab=\\,
cex=0.5,col=WD.col[WD.grp],
las=1,pch=16,
xaxs=\i\,yaxs=\i\
)
mtext(side=2,target.plot.var.title[k1],line=3)
mtext(side=1,\Air Temperature Day and Night\,line=2.8)
abline(h=0,col=\darkgrey\)
#axis(1, at = month.ticks, labels = FALSE, tcl = -0.3)
par(fig=c(0.7,1,0,1),new=T)
hist0<-hist(cdata[,target.plot.var[k1]],
plot=F,nclass=50)
hist1<-hist(cdata[,target.plot.var[k1]][WD.grp==1],
plot=F,breaks=hist0$breaks)
hist2<-hist(cdata[,target.plot.var[k1]][WD.grp==2],
plot=F,breaks=hist0$breaks)
barplot(hist1$counts,
axes=F,
horiz=T,
ylim=c(0,length(hist1$breaks)+1),
xlim=c(-5,max(c(hist1$counts,hist2$counts))),
space=0,col=WD.col[1],border=NA) # barplot
barplot(hist2$counts,
axes=F,
add=T,
horiz=T,
ylim=c(0,length(hist1$breaks)+1),
xlim=c(-5,max(c(hist1$counts,hist2$counts))),
space=0,col=WD.col[2],border=NA) # barplot
legend(0,
length(hist1$breaks)+1,
fill=WD.col,border=NA,
legend=WD.legend,bty=\n\,
cex=0.9)
dev.off()
}
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#water fluxes and cumulative sum of co2 and h20 fluxes
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-plot-begin eyJoZWlnaHQiOjMyOC44MDEsIndpZHRoIjo1MzIsInNpemVfYmVoYXZpb3IiOjAsImNvbmRpdGlvbnMiOltdfQ== -->
<img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAArwAAAGwCAMAAAB8TkaXAAAAZlBMVEUAAAAAADoAAGYAOpAAZrY6AAA6ADo6AGY6Ojo6OpA6kNtmAABmADpmAGZmOjpmOmZmtrZmtv+QOgCQZgCQ27aQ2/+2ZgC225C2/7a2///bkDrb////AAD/tmb/25D//7b//9v///+0e5yYAAAACXBIWXMAAA7DAAAOwwHHb6hkAAAT7klEQVR4nO2djXqjuBVASeLstPFsu0kbd0NCHL//SxYQYPwDugIJdMU5327imVwkHM7IF0lI2QlAKdnaJwAwFeQFtSAvqAV5QS3IC2pBXlAL8oJakBfUgrygFuQFtSAvqAV5QS3IC2pBXlAL8oJakBfUgrygFuQFtSAvqAV5QS3IC2pBXlAL8oJakBfUgrygFuQFtSAvqAV5QS3IC2pBXlAL8oJakBfUgrygFuQFtSAvqAV5QS3IC2pBXlAL8oJaPMubAcxmLXn9FgdbBHlBLcgLakFeUAvyglqQF9SCvKAW5AW1IC+oBXlBLcgLakFeUAvyglqQF9SCvKAW5AW1IC+oBXlBLVHK+1XitzpIkUDy/ryZh4weP5yL++pwqhG2Rxh58+zFvCjaF8LiOm2xF6wEkffnrVM2f/p0KO5CWfSFcYLIe9y/ti+Ly8TB6all7IVR4mp5ryH5hRFC5bxN0+ua816DvTBMoN6G494kBwPtrltx6At3ibKf9xqaX7iHCnlpfOEeSuQ94S/coEte7IUeeuQ90fjCJarkPSEv9NAlb6Wv19MAzeiW94tUYstok/ds79cdPJ0c6EChvMbRnrCD9qJz2qiT99T5evdHX1d/xN6E0Sfv2Gjxxd9Xr9E3ZbTKO/6zcxKBvQmjUd5RbhJg7E2W5OS9bWxpfFMlQXlvwN5E2YK86Jso25CXxDdJNiIvjW+KIC+oZTPynpiRlhzbkhd9k2JL8nLblhibkpfMIS02Jy/2psPG5KXpTQlXeYtmkcfX8Wh/9foGe9PBSd7jPtuZP/28DS967rde3yBvOrjIe/yz7+vln4LV6x3sTYat5bzImxDIC2pBXlAL8oJatigv+ibC9uRllC0Z5PLu5vaOTas3BNibBnJ5i8ePdOTF3hRwSBuOf/6diLzYmwZOI2z/+pWIvCQOSeB0w1bMm88wqd5AYG8CyOU9/tubuS71BgJ5E8BBXn8Jr0u9gUDeBNiqvNibABuWF3u1s1l5sVc/gW7Yft7M00KD3RPry8uq/+oJM7chz17Mi6J9Mau4QCCvclzlPe4tTWrFz1unbP70Oa/esGCvZhzl/XkbkPGC4757uHhoXAN5YTaO8va0HEFTy4u9enFueQdy2Evydl2HqHPeCuxVjGvO+/3HuyS6TY0Hk4xY5CVxUIyzvM+CGzaf9QYHe9XinDbsFq43OCQOaglyw3bcV5luEfkgRQv2aiXIDVstb93PcC171uF0lmFBXqUEuWGr5G20jbyrzIC9OnFOGzLBDVsl7/dzLW/kgxQGEgedBJnboK3lxV6dBJK3apx3p/bWbWZxC4C8Ggm1Yk7p78P78ABbbPJir0Y2uNzTXZgeqZCQ8o71q8Umr5majsCqQN6OL/xVhpu8ZnT4+1k0tUybvBXYqwnH3YBMx1f73YJCeU9kv4pwkvfQOjvYA3aBTnnp9FWDi7yCByT817sGNL46cFol0v5omv961wB5deAm77nlTVreE4MWKnDLebuZ6Iek04YTja8KnORtZor1XoSvdy2wN37c+nmbuQr54JQF//WuB/rGjvOTFO18saXqXRHsjRwm5oyBvVGDvKNgb8wg7zhMdogY5LWAvPGCvHbwN1KQ1w7yRsrUR9+zbN4whSZ5T4wWx4lry2uGKYrsVTYrcn69cUDjGyMTF9rLnz7nzYpUJi/2xsjEhfaKx495syK1yUvmECETF9qrWt6NyYu90TE55525Uq9KedE3Lpy7yur+hrLRPcybjq5QXuyNDfp5XUDfqEBeJ1iXJCaQ1xHsjQfknQD2xgHyTgF7owB5J0HqEANO6zb82e8du/xTsHojBXsjwHGhvfbZy5+3mXtgapcXfSPAfYTNx4TIFOTF3tUh550B+q4L8s4BeVcFeeeBvSsSSF6zso6OjbNnQdu7ImHk7RYzG9yILRV5yXtXJIi8giXUk5GXyQ7rEWRldEFgOvLS+K6Gs7ymVR2Xd1Mtbw32rkEQecuct2l60895DbS9axBG3m5xksHH4xOTF3vXIJC8/urVAvYuD/L6AnsXJ5C8mxmk6IG8SxPqhm0zgxQ9sHdh3OTtlogcn827va6yGuRdmCAjbBsbpOjA3mVheNgj3LMtS5Dh4e0NUjRg76IsPEhxTpqnnGz8YO+S0M/rFybpLAjyegZ5lwN5vYO9S4G83qHtXYog8gpGMxKWF3uXIkzL+/Nma5lTlpfEYSGCDA+ft7yaX69KsHcRQq3bUFgWhEpdXvRdgEAtr796dYK9SzCt5T0Mjfr6r1cr6BueKfKWLbBokcje/Jzp9aoFeYMzQd5CmjNsW15u24LjLm8+/EjwFci79hkkjqu8P2/ydHfr8mJvYBzllaa7NchL3hsUN3nF6a6/elWDvUFxkjfPZu30Pqle5WBvQBikCAz2hoNl/UODvcFA3uCQ+IYCecODvYFA3iXA3iAg7yJgbwiQdxnYdSUAzvIe99nT52HujMjNyYu9AXCVt3h4z58+j3vm87qDvp5xlLearlCtnZczSDEB5PWL88Sc11pe1m2YBvb6ZGLLe5DO6Z1db1qQOfhkWs6bO0yMnFlvYmCvR6b0NmTZw/ti9SYH9nqDft7FwV5fON+weZrRu2F5sdcXri3vwUvSsG156XPwxKRH30vmNsCblhd7/TAt57WvAumt3iQhcfACLe8qYK8PnAcpyHm9gLweoLdhJbB3PvTzrgRDbfNB3rVA3tm4rdvw0i3dQG/DfLB3JrS860HbOxPXG7bfpqfBNp+37pUYa6CR94S9c5kor+VJirxdB3Vju767gr2zcJL3cF6rbPQZtt7apvnArHXkrcHeOUxsecepdxs0DOUXyGvA3hkEuWGj5ZWDvdNxlbftK7PlvE3TS85rBXsn4zyf9+kz352+ny3PsLWODz6nibwd2DsV57kNL6eiWreBp4f9gbwTmbBuw/evj/r/ZerdAtg7jQnrNlQ9DuPymtWgCgYphJA4TMM1561GJw4vlrShlrcO6XWaTax3E2DvJJy7yg676m5svLOhkrfRlq4yEcg7hSD9vJW8TYcEgxQysHcCweSl5XWCuekTCLIPmwncndpbtzn1bgXsdSfUfN7S34f32wG2s/1uxW0B5HXFuZ+XJymCgb2OTGt5Z68tjbx3IHNwZGLacGDRkQBgrxsT5WVZ/zBgrwsT5R1/DEjQLYG898FeB6bJa1s4x74QH/IOgL1yJvY22GZE/rxZkmLkHQJ7xYTq5y0sW64g7yDYK4VFR+IDe4WEeYbNY71bBHtlOE9Gd3j+p/cQ8fR6Nwn2ipjwGJAY5J0M9kqY8BiQGOSdDvYKcM15v/+Qr+mPvDPAXjvO8j7Lb9iQdw7Ya8U5bWBPiqXAXhshb9i81LthmGRmIeQNm5d6t8wX+o4S8obNS73bBnvH4DGgyMHeYZjbEDvYOwgtb/Rg7xA8gBk/yDsAD2AqAHvvwwOYCqDP4T5BHsD0WS+csHeAMA9geqwXatD3DoEewPRXLxiQ9xb6edWAvdcgrxrIHK5xk9fMiPx+nj87B3ndwd4rnOQ97k2u235fol7ogb4XOMl7aJ0dXPDcf73QB3v7uMgr2A/bf71wCfqecduTonuMghG21UDeDjd5zy0v8q4H9ja45bzduNqBtGE9yBwanORtNgbsvQhfL9yyKXtH3qlbP2+zNVV+vUOVO8g7i+3oO/ZGnZ8ebrcHnAnyzmML8n71uRfA8LBWkm98W2uRNz1SX9RB8O6QVzEp2yt5b4HkNbnxyEPGyOuDdBtf0RsLI2/XHXGzc/ak4mCIZO0Vva0g8gomQSCvJ9JKfcd7F64JIq9gEgTy+kJ+rePH0jV2DS1vAqTir+ubCJXzNk0vOe8yJGGv81sI1NtgfcoYeX2jXl/386efNxmUJw8Tzh15E0KtvV/TBlwYpEgLjfY69C9cwiBFYqiz16Vv7Aq6ylJDl72zMh0GKZJDlb2zzpWWNz3U3LbNPVEGKRIkenvN6c3uHFl4kCLrcCoOHIlQ3+qEzs566ZSmnzdNYrP36x5zC0XeVIlKX6/OdgSSNy8TgzrtHVpaB3mDE5O9Yc4k0A3bw3uzbwXyrkg8+iqS13SV/byVt2vIuyaxzNUJdA5BBykOT5/Iuy5x2KtJ3m6Q4rBD3tVZTd+23i9V8nbJwnE/NK8MeZdjcXuvu3ND1R+st8EkDj9vyBsBy+p7Fjasu/TzboMF7Q0sbJ+Q8vbm5/goDuawlFBLNvLIuxUWsnfJ/AR5t0Rwe5e9NUTeTRG49V24XwN5t0VQe5fukqO3YXME83fxwRDk3R6hurIWH8hD3k0SojN2+UFo5N0qvv1dYQYF8m4aj/6uMPsHebdONw9hbimezscB5AUPGcQ60y6RF07z9V1nxjDyQsdUgdea7Y680GPS7PHVHtVAXrjC2d7VHpJDXriHg7/rPeCJvHAXefuLvBAfIn/XyncrkBdGsN2+LfRs0QDIC+OM6buqusgLEoYEXnkxHuQFCXfzh7UXkkJekPJ1wWnD8gLMZiV551YRKjaKkyDWPXaZgvxUgbzErlCQnyqQl9gVCvJTBfISu0JBfqpAXmJXKMhPFchL7AoF+akCeYldoSA/VSAvsSsU5KcK5CV2hYIAlgZ5QS3IC2pBXlAL8oJakBfUgrygFuQFtSAvqAV5QS3IC2pBXlAL8oJaQstbZNnDuyAub+MEB3w/Z9lOWri44O9fH+XXn7csy15s8Sb2uK+e096JYqvzyF5Hyz1X3R0kij0dnj4lsUXzYPnr2Pne/roEsYLfcRfidPEsBJa3KM+vEJxj/vhh4gQHFOW1Ou53ssIP5c+/n1/twcd9eQrlhS5D8sxSuIk9ff/R/FgQm1ch2dh5nKvuDhLFlsVW8spibb+42+sgiBVcvHOIy8WzEVZes8f2YWeLO+5fquCd5AATUv42JIXXBZ/yp09bcNkOVL7UntsKb2JPRf11/F02sdV7q0OGY7uqewcJYqv2vxRCFmu+D8feXgdBrODiXYUIL56VsPJe/uIsVG9OcEDX3kkKN/+2y6+W4CJ7Kc4/G4/vYvOd9UTa2E5e20nXDVJzkCS2+pf5n1JeWWwZ9mIt9+I6CGJvDxop1+XiWQksb32uhewUq09WwQHF49/7OoeTFN7Im73ag3s/O5QXbTTe/P3hH5ngRIqLtMF2HgfzEyOvJLYMqnJeYbkHQbkX10EQe3vQcKzTxbMSVt6u4ROEmpsK+wF5Vn/y7USFN//As1d78PkXWbZ8lsLr2OO+yjUPotju/sRyHkVzF1YfJImtPn4reWXlmixqNPbqOghibw8ajnW6eFaikbe6Dk+fEnkfmn+zosLNDZuTvEV7v2YV8iQ5kaaVLr9Wn9mW2Pa+SiCviS3TeYG8bbn9qMHfxcV1EMTeHjQc63bxbESUNtSppv0AkyeVTaqs8EN58/P3b0HB7c9MIyVJBWqsJ3KZvo7GFl3vlz1t6J2mNW3oyjWdarbfRf86CGJvDxqOdbx4FiK6Yavfk/2Aou2mkhf+/ct249GVW36wvdjPvC+v7UQuG9Gx2Lxz98r44di87bsVlWuyButV6V8HQeztQcOxEy7eCHF0lZm3Uoh6UI57eWyDoKusu7EyIwmWM++5ZT0RcWxbdXeQMLZuT0WxjWm2brX+aQpiBe+tC5lw8UaIZJCi+/0LDsgdYrv+cHtwI1n3wT0W3+SxdfeX7UQuc97h2HPV3UGyWJMMSGLbnw+f7+11EMQKLl4X4nLxrIQeHu6GBS0c2tFTwQFFOzQqiK2GcE2MLbj2pfkYrseNRuKbtEF00r3Y0ZPuVd3LYQSxTSYriO0+pYfP9/YtCWIFv4cuxOXi2WBiDqgFeUEtyAtqQV5QC/KCWpAX1IK8oBbkBbUgL6gFeUEtyAtqQV5QC/KCWpAX1IK8oBbkBbUgL6gFeUEtyAtqQV5QC/KCWpAX1IK8oBbkBbUgL6gFeUEtyAtqQV4hh24fKBFmOUQL3Y4Otwfc3VLrvCvW3FW+0gB5hTTrMj/PW5TzkmF5726p1Wxw1e3OtXmQV4iRd/aqnBcMy3tvS61mg6tudy5/56EV5BXSyNvsSWUSiNqpfFfv6ph1K4TXrysXj/u/9ubv82anzLy3b+Wp3cXstTa3/K8truNiJfV2gysvuzmkAfIKaeStv1eWmY2Zqq1tzHL2hbG3fW3krT77m/+rBrPe76m3LvSVvE1xvTqv9oToyyvOvhMGeYW08paN4PH3u7G0Eq+0q78tSPvayGt27Gu2gmg3kjpHX8nbFHcu63pLrf5OFciLvGJ68lbfijZvKJPP5kaqpn1t5DVetndXRsHevdaVvE1x7U9vt9QqejdsyIu8YvppQ5m6Pv6vzg8ePw4vTb9WI1Pz+kLepjXtbbvexl7K2xR3aqJvttQ67xBQ7c61zPuOGeQV0rthq9vO+svx939biQ69O6jy9XDL25A3W0H25e0Vd3dLrf4OWr9mbQKVBsgrpNdVZu6bqgb05+0f7cd8v6urdbHz8uUmpC6o6gc7J8av/eLubqnVk5eushPyimkHKV7ae7W6Zcy7PXfbtLR53Ze37mWoWtl689Kmia5KqgqqNjUtU41a3vN2f/e31Gq2PN6d6CmrQV4h/eHheny2lrBJSc/bSrWvL+S96Oc9b3FsSqs22/rLNNXnvob7W2qZlve8O9fWQd5ZfP/T66e35+JSB3lnkb/YY9YrLnWQdwbfz15vmzwXlz7IC2pBXlAL8oJakBfUgrygFuQFtSAvqAV5QS3IC2pBXlAL8oJakBfUgrygFuQFtSAvqAV5QS3IC2pBXlAL8oJakBfUgryglv8DiMQI+NnrM2cAAAAASUVORK5CYII=" />
<!-- rnb-plot-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#Work on (co)spectra data
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuaGVhZChjZGF0YSlcbmBgYFxuYGBgIn0= -->
```r
```r
head(cdata)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuXG5cbnpMLmN1dC5ybmc8LWMoLTUsLTAuMSwwLjEsMSkgIyMgc2V0dGluZyB0aGUgMyBaL0wgcmFuZ2VzIGZvciBwbG90dGluZ1xuY3NwZWMuZGF0YS5wcmU8LU5VTExcblxuZGF0YS5wcmU8LWNkYXRhICBcblxuIyMgd29yayBvbiBmaWxlIG5hbWUgb2Ygc3BlY3RydW0gZmlsZXNcbmNzcGVjLm5hbWUubHMucGFyc2U8LXBhc3RlKGNkYXRhJFRJTUVTVEFNUCR5ZWFyKzE5MDAsXG4gICAgICAgICAgICAgICAgICAgICAgICAgICBzcHJpbnRmKFxcJTAyZFxcLGNkYXRhJFRJTUVTVEFNUCRtb24rMSksXG4gICAgICAgICAgICAgICAgICAgICAgICAgICBzcHJpbnRmKFxcJTAyZFxcLGNkYXRhJFRJTUVTVEFNUCRtZGF5KSxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgIFxcLVxcLFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgc3ByaW50ZihcXCUwMmRcXCxjZGF0YSRUSU1FU1RBTVAkaG91ciksXG4gICAgICAgICAgICAgICAgICAgICAgICAgICBzcHJpbnRmKFxcJTAyZFxcLGNkYXRhJFRJTUVTVEFNUCRtaW4pLFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgXFxfYmlubmVkX2Nvc3BlY3RyYVxcLFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgc2VwPVxcXFwpXG5gYGBcbmBgYCJ9 -->
```r
```r
zL.cut.rng<-c(-5,-0.1,0.1,1) ## setting the 3 Z/L ranges for plotting
cspec.data.pre<-NULL
data.pre<-cdata
## work on file name of spectrum files
cspec.name.ls.parse<-paste(cdata$TIMESTAMP$year+1900,
sprintf(\%02d\,cdata$TIMESTAMP$mon+1),
sprintf(\%02d\,cdata$TIMESTAMP$mday),
\-\,
sprintf(\%02d\,cdata$TIMESTAMP$hour),
sprintf(\%02d\,cdata$TIMESTAMP$min),
\_binned_cospectra\,
sep=\\)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuIyMgc2NhbiBpZiBzcGVjdHJ1bSBmaWxlcyBleGlzdFxuY3NwZWMubmFtZS5sczwtcmVwKE5BLGxlbmd0aChjc3BlYy5uYW1lLmxzLnBhcnNlKSlcbmNzcGVjLmZpbGUuZXhpc3Q8LXJlcChGQUxTRSxsZW5ndGgoY3NwZWMubmFtZS5scy5wYXJzZSkpXG5jc3BlYy5maWxlLmxzPC1saXN0LmZpbGVzKHBhc3RlMChlZGR5cHJvLnBhdGgsXFxlZGR5cHJvX2Jpbm5lZF9jb3NwZWN0cmFcXFxcXFwpKVxuZm9yKGkgaW4gMTpsZW5ndGgoY3NwZWMuZmlsZS5leGlzdCkpe1xuICBjc3BlYy5uYW1lLmxzW2ldPC1jc3BlYy5maWxlLmxzW3doaWNoKGdyZXBsKGNzcGVjLm5hbWUubHMucGFyc2VbaV0sXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY3NwZWMuZmlsZS5scykpWzFdXVxuICBjc3BlYy5maWxlLmV4aXN0W2ldPC1pZmVsc2UoIWlzLm5hKGNzcGVjLm5hbWUubHNbaV0pLFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgVFJVRSxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIEZBTFNFKVxufVxuXG5jc3BlYy52YXIubmFtZTwtYyhcXG4uZlxcLFxcZlxcLFxcZmgudVxcLFxcZlN1XFwsXFxmU3ZcXCxcXGZTd1xcLFxcZlN0XFwsXFxmU2NcXCxcXGZTcVxcLFxcZlNtXFwsXFxmU25cXCxcbiAgICAgICAgICAgICAgICAgIFxcZkN3dVxcLFxcZkN3dlxcLFxcZkN3dFxcLFxcZkN3Y1xcLFxcZkN3cVxcLFxcZkN3bVxcLFxcZkN3blxcKVxuY3NwZWMudmFyLm5hbWUuc2VsZWN0PC1jKFxcbi5mXFwsXFxmXFwsXFxmaC51XFwsXFxmU3VcXCxcXGZTdlxcLFxcZlN3XFwsXFxmU3RcXCxcXGZTY1xcLFxcZlNxXFwsXG4gICAgICAgICAgICAgICAgICAgICAgICAgXFxmQ3d1XFwsXFxmQ3d2XFwsXFxmQ3d0XFwsXFxmQ3djXFwsXFxmQ3dxXFwpXG5cbmBgYFxuYGBgIn0= -->
```r
```r
## scan if spectrum files exist
cspec.name.ls<-rep(NA,length(cspec.name.ls.parse))
cspec.file.exist<-rep(FALSE,length(cspec.name.ls.parse))
cspec.file.ls<-list.files(paste0(eddypro.path,\eddypro_binned_cospectra\\\))
for(i in 1:length(cspec.file.exist)){
cspec.name.ls[i]<-cspec.file.ls[which(grepl(cspec.name.ls.parse[i],
cspec.file.ls))[1]]
cspec.file.exist[i]<-ifelse(!is.na(cspec.name.ls[i]),
TRUE,
FALSE)
}
cspec.var.name<-c(\n.f\,\f\,\fh.u\,\fSu\,\fSv\,\fSw\,\fSt\,\fSc\,\fSq\,\fSm\,\fSn\,
\fCwu\,\fCwv\,\fCwt\,\fCwc\,\fCwq\,\fCwm\,\fCwn\)
cspec.var.name.select<-c(\n.f\,\f\,\fh.u\,\fSu\,\fSv\,\fSw\,\fSt\,\fSc\,\fSq\,
\fCwu\,\fCwv\,\fCwt\,\fCwc\,\fCwq\)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuIyMgd2hpY2ggc3BlY3RydW0gdmFyaWFibGVzIHRvIHBsb3QgXG50YXJnZXQuY3NwZWM8LWMoXFxmU3dcXCxcXGZTdFxcLFxcZlNjXFwsXFxmU3FcXCxcbiAgICAgICAgICAgICAgICBcXGZDd3VcXCxcXGZDd3RcXCxcXGZDd2NcXCxcXGZDd3FcXClcbmBgYFxuYGBgIn0= -->
```r
```r
## which spectrum variables to plot
target.cspec<-c(\fSw\,\fSt\,\fSc\,\fSq\,
\fCwu\,\fCwt\,\fCwc\,\fCwq\)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
#for running just one month comment out 575 and 659
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuIyMjIGxvb3AgdGhyb3VnaCB0aGUgc3BlY3RydW0gZmlsZXMsIGNvbWJpbmUgdGhlbSBpbnRvIGEgc2luZ2xlIGRhdGFmcmFtZSBmb3IgcGxvdHRpbmdcbm1vbnRoLmlkPC1wYXN0ZTAoY2RhdGEkVElNRVNUQU1QJHllYXIrMTkwMCxcXC1cXCxzcHJpbnRmKFxcJTAyZFxcLGNkYXRhJFRJTUVTVEFNUCRtb24rMSkpXG5tb250aC5pZC5sczwtbGlzdChuYW1lcyh0YWJsZShtb250aC5pZCkpLFxuICAgICAgICAgICAgICAgICAgYXMudmVjdG9yKHRhYmxlKG1vbnRoLmlkKSkpXG4jY2FuIGNvbW1lbnQgb3V0IDU3NSBhbmQgNjU5IGlmIHlvdSB3YW50IHRvIHJ1biBqdXN0IG9uZSBtb250aFxuZm9yKGoxIGluIDE6bGVuZ3RoKG1vbnRoLmlkLmxzW1sxXV0pKXtcbiNqMT03ICNhZGQganVzdCB0aGUgbW9udGggaGVlcmUuIGNvbW1lbnQgb3V0IGlmIHdhbiB0byBydW4gZnVsbCBsb29wXG5jc3BlYy5kYXRhLnByZTwtZGF0YS5mcmFtZSgpXG50YXJnZXQubW9uPC1tb250aC5pZC5sc1tbMV1dW2oxXVxudGFyZ2V0Lm1vbi5sb2M8LXdoaWNoKG1vbnRoLmlkPT10YXJnZXQubW9uKVxuXG5wcmludChwYXN0ZTAoXFwjIyMjIyMgIFxcLHRhcmdldC5tb24sXFwgICMjIyMjI1xcKSlcblxuZm9yKGkxIGluIDE6bGVuZ3RoKGNzcGVjLm5hbWUubHNbdGFyZ2V0Lm1vbi5sb2NdKSl7XG4gIGlmKGNzcGVjLmZpbGUuZXhpc3RbdGFyZ2V0Lm1vbi5sb2NdW2kxXSl7XG4gICAgXG4gICAgIyMgcmVhZCBpbiBlYWNoIChjbylzcGVjdHJhIGZpbGVzXG4gICAgY3NwZWMuZGF0YTwtcmVhZC5jc3YocGFzdGUoZWRkeXByby5wYXRoLFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIFxcZWRkeXByb19iaW5uZWRfY29zcGVjdHJhXFxcXFxcLFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGNzcGVjLm5hbWUubHNbdGFyZ2V0Lm1vbi5sb2NdW2kxXSxzZXA9XFxcXCksXG4gICAgICAgICAgICAgICAgICAgICAgICAgc2tpcD0xMSxcbiAgICAgICAgICAgICAgICAgICAgICAgICBoZWFkZXI9VCxcbiAgICAgICAgICAgICAgICAgICAgICAgICBuYS5zdHJpbmdzPWMoXFwtOTk5OVxcLFxcLTk5OTkuMFxcKSlcbiAgICBcbiAgICBjb2xuYW1lcyhjc3BlYy5kYXRhKTwtY3NwZWMudmFyLm5hbWVcbiAgICBcbiAgICAjIyBjb21iaW5lIHNwZXRydW0gZGF0YSB3aXRoIGJhc2ljIHN0YXRlIHZhcmlhYmxlcywgZS5nLiwgV1MsIFdEXG4gICAgY3NwZWMuZGF0YS50bXA8LWRhdGEuZnJhbWUoZGF0ZT1yZXAoZGF0YS5wcmUkZGF0ZVt0YXJnZXQubW9uLmxvY11baTFdLG5yb3coY3NwZWMuZGF0YSkpLFxuICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpbWU9cmVwKGRhdGEucHJlJHRpbWVbdGFyZ2V0Lm1vbi5sb2NdW2kxXSxucm93KGNzcGVjLmRhdGEpKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB6TD1yZXAoZGF0YS5wcmUkWC56LmQuLkxbdGFyZ2V0Lm1vbi5sb2NdW2kxXSxucm93KGNzcGVjLmRhdGEpKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBMPXJlcChkYXRhLnByZSRMW3RhcmdldC5tb24ubG9jXVtpMV0sbnJvdyhjc3BlYy5kYXRhKSksXG4gICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgV1M9cmVwKGRhdGEucHJlJHdpbmRfc3BlZWRbdGFyZ2V0Lm1vbi5sb2NdW2kxXSxucm93KGNzcGVjLmRhdGEpKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBXRD1yZXAoZGF0YS5wcmUkd2luZF9kaXJbdGFyZ2V0Lm1vbi5sb2NdW2kxXSxucm93KGNzcGVjLmRhdGEpKSxcbiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBjc3BlYy5kYXRhWyxjc3BlYy52YXIubmFtZS5zZWxlY3RdKVxuICAgIFxuICAgICMjIyBmaWx0ZXIgc3BldHJ1bSBkYXRhIGJhc2VkIG9uIGNvcnJlc3BvbmRpbmcgdmFyaWFuY2UvY292YXJpYW5jZVxuICAgIFxuICAgIGlmKGlzLm5hKGRhdGEucHJlJHUuW3RhcmdldC5tb24ubG9jXVtpMV0pKXtcbiAgICAgIGNzcGVjLmRhdGEudG1wJGZTdTwtTkFcbiAgICAgIGNzcGVjLmRhdGEudG1wJGZDd3U8LU5BXG4gICAgICBjc3BlYy5kYXRhLnRtcCRmU3Y8LU5BXG4gICAgICBjc3BlYy5kYXRhLnRtcCRmQ3d2PC1OQVxuICAgICAgY3NwZWMuZGF0YS50bXAkZlN3PC1OQVxuICAgIH1cbiAgICBpZihpcy5uYShkYXRhLnByZSRIW3RhcmdldC5tb24ubG9jXVtpMV0pKXtcbiAgICAgIGNzcGVjLmRhdGEudG1wJGZTdDwtTkFcbiAgICAgIGNzcGVjLmRhdGEudG1wJGZDd3Q8LU5BXG4gICAgfVxuICAgIGlmKGlzLm5hKGRhdGEucHJlJGNvMl9mbHV4W3RhcmdldC5tb24ubG9jXVtpMV0pKXtcbiAgICAgIGNzcGVjLmRhdGEudG1wJGZTYzwtTkFcbiAgICAgIGNzcGVjLmRhdGEudG1wJGZDd2M8LU5BXG4gICAgfVxuICAgIGlmKGlzLm5hKGRhdGEucHJlJExFW3RhcmdldC5tb24ubG9jXVtpMV0pKXtcbiAgICAgIGNzcGVjLmRhdGEudG1wJGZTcTwtTkFcbiAgICAgIGNzcGVjLmRhdGEudG1wJGZDd3E8LU5BXG4gICAgfVxuICAgIGNzcGVjLmRhdGEucHJlPC1yYmluZC5kYXRhLmZyYW1lKGNzcGVjLmRhdGEucHJlLGNzcGVjLmRhdGEudG1wKVxuICB9ZWxzZXtcbiAgICBwcmludChwYXN0ZShjc3BlYy5uYW1lLmxzLnBhcnNlW3RhcmdldC5tb24ubG9jXVtpMV0sXFxoYXMgbm8gc3BlY3RydW0gZmlsZXNcXCkpICBcbiAgfSAgICAgICAgICAgICAgICAgICAgICAgXG59XG5cblxuIyMjIG91dHB1dCBzaXRlLXNwZWNpZmljIChjbylzcGVjdHJhIGZpbGUgKGNvbXBvc2l0ZSBhbGwgcmVjb3JkcylcbndyaXRlLmNzdihjc3BlYy5kYXRhLnByZSxcbiAgICAgICAgICBwYXN0ZTAocGxvdC5wYXRoLFxuICAgICAgICAgICAgICAgICB0YXJnZXQubW9uLFxuICAgICAgICAgICAgICAgICBcXF9jb3NwZWN0cmFfY29tcGlsZWRfXFwsXG4gICAgICAgICAgICAgICAgIFN5cy5EYXRlKCksXFwuY3N2XFwsXG4gICAgICAgICAgICAgICAgIHNlcD1cXFxcKSxcbiAgICAgICAgICByb3cubmFtZXM9RilcblxuIyMjIyBsb29wIHRocm91Z2ggdGhlIHRhcmdldC52YXIsIGRvIHNwZWN0cnVtIHBsb3R0aW5nIFxuZm9yKG0yIGluIDE6bGVuZ3RoKHRhcmdldC5jc3BlYykpe1xuICBcbiAgY29zcGVjdHJhX3Bsb3QzKGNzcGVjLmRhdGEucHJlPWNzcGVjLmRhdGEucHJlLFxuICAgICAgICAgICAgICAgICAgdGFyZ2V0LnZhcj10YXJnZXQuY3NwZWNbbTJdLCBcbiAgICAgICAgICAgICAgICAgIGNhc2U9XFxDb25jb3JkXFwsXG4gICAgICAgICAgICAgICAgICB5ZWFyPWNkYXRhJFRJTUVTVEFNUCR5ZWFyW3RhcmdldC5tb24ubG9jXVsxXSsxOTAwLFxuICAgICAgICAgICAgICAgICAgZG95Lmk9Y2RhdGEkVElNRVNUQU1QJHlkYXlbdGFyZ2V0Lm1vbi5sb2NdWzFdKzEsXG4gICAgICAgICAgICAgICAgICBkb3kuZj1jZGF0YSRUSU1FU1RBTVAkeWRheVt0YXJnZXQubW9uLmxvY11bbGVuZ3RoKHRhcmdldC5tb24ubG9jKV0rMSxcbiAgICAgICAgICAgICAgICAgIG91dHB1dD1ULFxuICAgICAgICAgICAgICAgICAgb3V0RGlyPXBsb3QucGF0aCxcbiAgICAgICAgICAgICAgICAgIHBsb3QubG9lc3M9VCxcbiAgICAgICAgICAgICAgICAgIHBvc3RmaXg9cGFzdGUwKFxcX1xcLFN5cy5EYXRlKCkpLFxuICAgICAgICAgICAgICAgICAgbG9nLnkudmFsdWU9VCxcbiAgICAgICAgICAgICAgICAgIHpMLmN1dC5ybmc9ekwuY3V0LnJuZykgIFxuICBcbn1cbn0gI2NvbW1lbnQgb3V0IGZvciBqdXN0IHJ1bm5pbmcgb25lIG1vbnRoLiByZW1vdmUgIyBpZiB5b3Ugd2FudCB0byBydW4gYSBmdWxsIGxvb3AuIFxuXG5cblxuYGBgXG5gYGAifQ== -->
```r
```r
### loop through the spectrum files, combine them into a single dataframe for plotting
month.id<-paste0(cdata$TIMESTAMP$year+1900,\-\,sprintf(\%02d\,cdata$TIMESTAMP$mon+1))
month.id.ls<-list(names(table(month.id)),
as.vector(table(month.id)))
#can comment out 575 and 659 if you want to run just one month
for(j1 in 1:length(month.id.ls[[1]])){
#j1=7 #add just the month heere. comment out if wan to run full loop
cspec.data.pre<-data.frame()
target.mon<-month.id.ls[[1]][j1]
target.mon.loc<-which(month.id==target.mon)
print(paste0(\###### \,target.mon,\ ######\))
for(i1 in 1:length(cspec.name.ls[target.mon.loc])){
if(cspec.file.exist[target.mon.loc][i1]){
## read in each (co)spectra files
cspec.data<-read.csv(paste(eddypro.path,
\eddypro_binned_cospectra\\\,
cspec.name.ls[target.mon.loc][i1],sep=\\),
skip=11,
header=T,
na.strings=c(\-9999\,\-9999.0\))
colnames(cspec.data)<-cspec.var.name
## combine spetrum data with basic state variables, e.g., WS, WD
cspec.data.tmp<-data.frame(date=rep(data.pre$date[target.mon.loc][i1],nrow(cspec.data)),
time=rep(data.pre$time[target.mon.loc][i1],nrow(cspec.data)),
zL=rep(data.pre$X.z.d..L[target.mon.loc][i1],nrow(cspec.data)),
L=rep(data.pre$L[target.mon.loc][i1],nrow(cspec.data)),
WS=rep(data.pre$wind_speed[target.mon.loc][i1],nrow(cspec.data)),
WD=rep(data.pre$wind_dir[target.mon.loc][i1],nrow(cspec.data)),
cspec.data[,cspec.var.name.select])
### filter spetrum data based on corresponding variance/covariance
if(is.na(data.pre$u.[target.mon.loc][i1])){
cspec.data.tmp$fSu<-NA
cspec.data.tmp$fCwu<-NA
cspec.data.tmp$fSv<-NA
cspec.data.tmp$fCwv<-NA
cspec.data.tmp$fSw<-NA
}
if(is.na(data.pre$H[target.mon.loc][i1])){
cspec.data.tmp$fSt<-NA
cspec.data.tmp$fCwt<-NA
}
if(is.na(data.pre$co2_flux[target.mon.loc][i1])){
cspec.data.tmp$fSc<-NA
cspec.data.tmp$fCwc<-NA
}
if(is.na(data.pre$LE[target.mon.loc][i1])){
cspec.data.tmp$fSq<-NA
cspec.data.tmp$fCwq<-NA
}
cspec.data.pre<-rbind.data.frame(cspec.data.pre,cspec.data.tmp)
}else{
print(paste(cspec.name.ls.parse[target.mon.loc][i1],\has no spectrum files\))
}
}
### output site-specific (co)spectra file (composite all records)
write.csv(cspec.data.pre,
paste0(plot.path,
target.mon,
\_cospectra_compiled_\,
Sys.Date(),\.csv\,
sep=\\),
row.names=F)
#### loop through the target.var, do spectrum plotting
for(m2 in 1:length(target.cspec)){
cospectra_plot3(cspec.data.pre=cspec.data.pre,
target.var=target.cspec[m2],
case=\Concord\,
year=cdata$TIMESTAMP$year[target.mon.loc][1]+1900,
doy.i=cdata$TIMESTAMP$yday[target.mon.loc][1]+1,
doy.f=cdata$TIMESTAMP$yday[target.mon.loc][length(target.mon.loc)]+1,
output=T,
outDir=plot.path,
plot.loess=T,
postfix=paste0(\_\,Sys.Date()),
log.y.value=T,
zL.cut.rng=zL.cut.rng)
}
} #comment out for just running one month. remove # if you want to run a full loop.
```